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CHAPTER  I 


INTRODUCTION 


As  the  demand  for  more  precise  control  of  material  parameters  such 
as  thermal  expansion,  directional  strength  and  stiffness,  and  impact  and 
fatigue  resistance  in  structures  has  increased,  designers  have  begun  to 
move  away  from  metal  components  and  are  now  attempting  to  "tailor"  com¬ 
posite  materials  to  suit  their  particular  needs.  Thus,  the  use  of  fiber- 
reinforced  composites  such  as  graphite/epoxy  and  glass/epoxy  has  increased 
rapidly  over  the  past  ten  years. 

While  an  extensive  effort  during  the  past  decade  has  been  devoted 
to  studying  the  above  topics,  it  has  only  been  recently  that  viscoelastic, 
i.e.,  time-dependent,  behavior  of  composites  has  been  investigated.  The 
viscoelastic  behavior  of  a  material  can  be  of  paramount  importance  when 
a  composite  component  is  subjected  to  repeated  loading  cycles  or  sub¬ 
stantial  loads  at  moderately  high  temperatures  for  long  periods  of  time. 
Environmental  cycling  or  fluctuations  of  the  ambient  temperature  and/or 
relative  humidity  can  also  have  considerable  influence  on  the  performance 
of  the  composite.  Under  such  conditions,  the  internal  stress  distribu¬ 
tions  and  overall  strain  of  the  composite  could  ultimately  change  enough 
with  time  to  cause  failure  of  the  component.  Alternatively,  the  stress 
state  could  also  change  so  as  to  relieve  high  stress  states  in  certain 


areas  of  the  composite,  preventing  failure  of  the  component. 

It  is  not  possible  to  accurately  simulate  environmental  or  loading 


2 


cycles  with  a  time-independent  analysis.  This  approach  does  not  realisti¬ 
cally  model  the  situation  since  it  neglects  any  redistribution  or  changes 
in  the  internal  stress  state  with  time.  And,  since  the  distribution  of 
stresses  can  change  one  way  or  another  with  time,  the  ability  to  predict 
these  changes  can  be  extremely  important . 

As  an  example,  graphite/epoxy  or  other  composite  components 
on  aircraft  are  often  exposed  to  significant  environmental  changes. 

Whether  these  are  detrimental  to  the  performance  of  the  composite  or  not 
can  be  determined  by  means  of  a  time-dependent  analysis.  Spacecraft  and 
satellite  components  must  be  dimensionally  stable  throughout  their  service 
lifetime,  to  preserve  the  functional  integrity  of  the  optical  or  electronic 
hardware  which  they  carry.  During  fabrication,  storage,  and  launch,  the 
satellite  sees  a  variety  of  environmental  fluctuations,  which  cause  di¬ 
mensional  changes  to  the  structure  and  the  Instruments  it  carries.  To 
accurately  determine  these  dimensional  changes,  a  time-dependent  analysis 
should  be  used.  Automotive  components  such  as  driveshafts  and  engine 
motor  mounts  must  be  able  to  withstand  sustained  loading  and/or  cyclic 
loadings  for  long  periods  of  time.  Whether  the  composite  is  capable  of 
such  performance  may  also  be  determined  by  means  of  a  time-dependent 
analysis. 

In  order  to  experimentally  characterize  the  time-dependent  behavior 
of  a  composite  material,  a  number  of  tests  must  be  conducted.  There  are 
many  factors  (e.g.,  temperature,  moisture,  loading  rate,  and  stress  level) 
that  affect  the  viscoelastic  properties  of  a  composite.  Thus,  an  exten¬ 
sive  amount  of  testing  must  be  undertaken  to  determine  how  and  to  what 
extent  each  of  these  factors  will  affect  the  composite.  While  a  particu- 
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lar  composite  material  system  might  be  adequately  characterized  by  such 
a  testing  program,  If  a  designer  wishes  to  use  a  different  fiber-matrix 
combination,  he  must  re-evaluate  all  the  properties  determined  in  the 
previous  testing  program.  This  can  be  a  problem  since  there  are  numer¬ 
ous  fiber-matrix  combinations  available  to  a  designer. 

To  overcome  this  difficulty,  a  micromechanics  analysis  and  asso¬ 
ciated  computer  code  has  been  developed  in  the  present  study  which  in¬ 
cludes  time-dependent  effects.  By  inputing  the  appropriate  constituent 
viscoelastic  properties  to  the  program,  the  designer  can  try  any  combin¬ 
ation  of  fiber  and  matrix  material,  to  determine  if  that  composite  will 
satisfy  his  need.  Since  only  the  individual  constituent  properties  need 
to  be  determined  experimentally,  this  analysis  becomes  extremely  cost 
effective  considering  the  amount  of  time  and  experimental  testing  saved 
by  not  having  to  evaluate  all  fiber-matrix  combinations.  Also,  most 
fibers  do  not  exhibit  time-dependence.  Therefore,  only  the  candidate 
matrix  materials  need  be  tested. 

To  model  the  various  fiber-matrix  combinations,  the  analysis  uses 
a  finite  element  model  capable  of  simulating  a  unidirectional  composite 
subjected  to  any  combination  of  longitudinal  and  transverse  normal  load¬ 
ings,  as  well  as  hygrothermal  loading.  Time-independent  nonlinear 
(elastoplastic)  material  behavior  is  included,  as  is  a  creep  formulation 
which  uses  stresses  as  the  independent  state  variable.  The  program  also 
contains  two  failure  criteria,  viz.,  an  octahedral  shear  stress  criterion 
and  a  hydrostatic  criterion. 

Since  the  effects  of  temperature  and  moisture  are  included  in  the 
analysis,  it  is  possible  to  examine  the  internal  stress  state  of  a  com- 
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posite  as  a  function  of  time  after  curing  a  laminate,  as  well  as  during 
external  loading.  The  analysis  can  also  model  any  arbitrary  multiple 
step  loading  (stress  reversal,  cycling,  etc.)  as  well  as  recovery. 

Some  numerical  examples  of  the  capability  are  presented  in  Chapter  V. 


CHAPTER  II 


LITERATURE  SURVEY 

Whether  considering  creep  or  relaxation  there  are  numerous  theories 
available,  each  having  its  strong  and  weak  points.  A  review  of  various 
theories  will  be  presented  in  this  chapter  and  the  relative  advantages 
and  disadvantages  of  each  will  be  discussed.  Before  this  discussion  be¬ 
gins,  a  few  of  the  terms,  viz.,  creep,  recovery,  and  relaxation  will  be 
briefly  defined. 

In  a  creep  test,  stress  is  the  independent  state  variable,  while  in 
a  relaxation  test,  strain  is  the  independent  state  variable.  Recovery 
occurs  when  all  of  the  stress  is  removed  after  a  creep  test.  There  are 
usually  three  stages  of  creep:  primary,  during  which  a  rapid  decrease 
of  the  creep  rate  occurs,  secondary,  where  the  creep  rate  is  more  or  less 
constant,  and  tertiary,  where  a  rapid  increase  in  creep  rate  occurs.  The 
tertiary  stage  of  creep  is  usually  brief,  leading  quickly  to  failure  and 
is  not  cf  much  use  in  design.  Illustrative  examples  of  creep-recovery 
and  relaxation  are  shown  in  Figures  la  and  lb,  respectively.  The  three 
stages  of  creep  are  shown  in  Figure  1c. 

Empirical  Approaches 

Observation  of  the  time-dependence  or  creep  of  metals  was  first 
reported  around  the  mid-nineteenth  century  [1].  Early  investigators 
used  an  empirical  approach  in  dealing  with  the  observed  time-dependence 
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a)  Plot  of  strain  versus  time  for  o  ■  constant 
(creep)  and  c  =  zero  (recovery) 
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b)  Plot  of  stress  versus  time  for  e  »  constant 
(relaxation) 


c)  Plot  of  strain  versus  time  showing  the  three 
stages  of  creep 

Illustrative  examples  of  creep-recovery,  relaxation, 
and  the  three  stages  of  creep 


of  metals  [1,2].  Many  equations  were  written  which  yielded  good  agree¬ 
ment  with  experimental  data  for  uniaxial  creep  at  low  stresses  under  a 
constant  externally  applied  load.  These  equations  were  found  to  be  ade¬ 
quate  to  describe  either  primary  or  secondary  stages  of  creep,  but  rarely 
both.  Difficulties  arose  when  investigators  tried  to  describe  combined 
states  of  stress  using  this  empirical  approach.  The  multiaxial  stress 
experiments  necessary  were  not  easily  performed,  and  the  mathematical 
representation  became  more  involved  than  those  for  uniaxial  creep.  Also, 
if  the  situation  involved  a  changing  stress  history  or  a  loading  that  was 
a  function  of  time,  the  equation  for  the  creep  rate  had  to  be  a  function 
of  strain  as  well  as  the  applied  stress  in  order  to  fit  the  experimental 
data.  This  fact,  along  with  the  difficulty  of  multiaxial  representation, 
caused  many  to  realize  that  a  more  rigorous  approach  was  necessary. 

Mechanical  Models 

A  number  of  investigators  have  utilized  mechanical  models  to  repre¬ 
sent  the  creep  behavior  of  materials  [1,3].  By  assembling  springs  and 
dashpots  in  a  variety  of  configurations,  they  have  been  able  to  describe 
many  different  types  of  materials;  from  stainless  steels,  which  exhibit 
a  great  deal  of  secondary  creep,  to  polymers,  which  do  not.  These  models 
become  extremely  difficult  to  visualize,  much  less  formulate,  when  a 
combined  stress  state  is  considered. 

Differential  Operator  Form 


Both  of  the  time-dependent  stress-strain  formulations  just  discus¬ 
sed  were  primarily  phenomenological,  models  being  constructed  and  equa- 
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tlons  written  to  describe  specific  experimental  data.  In  creep,  actual 
values  of  strain  (response)  depend  upon  the  entire  stress  (input)  his¬ 
tory  and  not  just  the  current  value  of  stress.  That  is,  the  creep  be¬ 
havior  of  most  materials,  composites  included,  are  affected  by  the 
total  loading  history.  In  order  to  account  for  this,  a  more  mathematical 
representation  is  needed. 

One  of  these  mathematical  forms  is  the  differential  operator  method 
[1,3].  A  differential  equation  is  written 

P a  =  Qe  (1) 

where  P  and  Q  are  a  series  of  linear  differential  operators  with  respect 
to  time.  By  choosing  the  proper  terms  in  P  and  Q  and  excluding  others, 
any  material  behavior  may  be  represented.  It  can  be  shown  that  Eq.  (1) 
can  be  reduced  to  give  a  representation  of  any  of  the  mechanical  model 
combinations.  A  generalization  to  multlaxial  stress  states  is  given  by 
Findley  and  Lai  [1]. 

A  common  solution  technique  for  Eq.  (1)  utilizes  the  Correspondence 
Principle  [1].  This  method  utilizes  a  known  linearly  elastic  (time-inde¬ 
pendent)  solution,  selected  in  accordance  with  the  particular  problem.  The 
Laplace  transform  of  this  solution  is  then  taken  to  obtain  equations  of 
state  in  the  imaginary  or  s-domain.  Next,  a  substitution  is  made  for 
the  modulus  and  Poisson's  ratio,  e.g.,  E  ■  sE(s)  and  v  "  sv(s),  so  that 
these  quantities  will  become  time-dependent  terms  when  the  inverse  Lap¬ 
lace  transform  is  taken.  Following  this  substitution,  the  desired  quan¬ 
tities  are  solved  for.  The  solution  in  the  real  (time)  domain  is  then 
generated  by  taking  the  inverse  Laplace  transform  of  the  resulting  equa¬ 
tion.  The  limitation  in  using  this  solution  technique  is  that  Eq.  (1)  is 
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tion.  The  limitation  in  using  this  solution  technique  is  that  Eq.  (1)  is 


now  restricted  to  the  linear  viscoelastic  range  since  linearly  elastic 
(time-independent)  equations  are  used.  This  restriction  means  that 
the  solution  is  only  applicable  to  a  material  which  is  subjected  to 
stresses  (or  strains)  that  are  sufficiently  low.  This  restriction  will 
be  discussed  in  more  detail  in  the  next  chapter.  Of  course,  it  is  possi 
ble  to  find  a  solution  without  using  the  Correspondence  Principle,  but 
the  differential  equations  can  become  cumbersome  and  extremely  difficult 
to  solve. 


Integral  Representation 

Another  and  more  widely  used  mathematical  model  is  the  integral 
representation  [1,3-8],  The  advantages  of  this  method  over  the  differ¬ 
ential  operator  form  are  the  relative  ease  of  incorporating  measured  ex¬ 
perimental  data  and  the  inclusion  of  temperature  and  moisture  effects. 
For  these  reasons  this  method  has  been  very  appealing  to  many  investi¬ 
gators,  and  consequently  is  becoming  increasingly  popular. 

Probably  the  simplest  of  the  integral  representations  is  the 
Boltzmann  superposition  principle  [1].  It  states  that  the  sum  of  the 
strain  (stress)  outputs  caused  by  several  stress  (strain)  inputs  equals 
the  sum  of  the  strain  (stress)  outputs  caused  by  each  stress  (strain) 
input  individually.  Thus,  this  representation  is  restricted  to  linear 
viscoelastic  behavior,  and  as  discussed  before,  can  only  be  used  when 
the  applied  stresses  (or  strains)  are  sufficiently  low. 

A  slightly  more  general  integral  form  was  proposed  by  Leaderman  [9] 


This  representation,  referred  to  as  the  modified  superposition  princi¬ 
ple,  modified  the  Botlzmann  superposition  integral  slightly  so  that 
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stress-dependence  could  be  accounted  for.  For  Instance,  the  creep  com¬ 
pliance  is  usually  divided  into  two  portions,  i.e.,  a  time-dependent 
portion  and  a  time-independent  portion.  While  this  method  allows  for 
the  stress-dependence  of  the  time-dependent  portion,  it  does  not  for  the 
time-independent  portion.  This  method  has  been  useful;  however,  it 
lacks  the  complete  generality  necessary  to  describe  nonlinear  behavior. 

Lai  and  Findley  [10]  have  recently  used  the  modified  superposition 
principle  in  conjunction  with  other  theories  to  represent  the  nonlinear 
viscoelastic  behavior  of  aluminum.  In  this  study,  multiaxial  creep 
tests  were  run  to  verify  their  hybrid  theory  with  experimental  data. 
Excellent-to-fair  agreement  was  obtained  between  the  measured  data  and 
the  predictions  based  on  constant-stress  tests.  This  relatively  new 
theory,  while  having  potential  applications,  has  not  been  applied  to 
many  different  types  of  materials.  Also,  it  is  somewhat  empirical  In 
nature  and  requires  several  experimental  tests  to  obtain  the  parameters 
used  in  the  constitutive  equations. 

Multiple  integral  representations  [5,6]  are  very  attractive  in 
that  they  are  completely  general  In  nature  and  can  account  for  any  de¬ 
gree  of  nonlinearity.  Hence,  they  are  not  limited  to  any  one  type  or 
class  of  materials.  The  difficulty  with  this  method  arises  when  the 
material  contains  a  high  degree  of  nonlinearity,  making  simplification 
of  many  terms  Impossible.  Without  such  simplification,  the  experimental 
requirements  of  this  theory  are  overwhelming  and  the  kernels,  i.e.,  the 
stress-, temperature-,  moisture-dependent  parameters  in  the  equations, 
extremely  difficult  to  evaluate  experimentally. 

Schapery  has  developed  much  simpler,  single-integral  nonlinear  con- 
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stitutive  equations,  derived  from  thermodynamic  theory  in  terms  of 
strains  [7]  or  stresses  [8].  These  equations  contain  kernels  which 
are  written  in  terms  of  linear  viscoelastic  properties,  and  are  much 
easier  to  evaluate  than  those  in  the  multiple  integral  forms.  In  fact, 
the  equations  can  be  reduced  to  the  Boltzmann  superposition  principle  by 
equating  certain  terms  to  one;  and  experimental  determination  of  the 
kernels  requires  only  uniaxial  creep  tests  followed  by  a  recovery  period. 

For  these  reasons,  the  theory  developed  by  Schapery  [8]  was  chosen 
for  use  in  the  present  analysis.  This  decision  was  based  upon  the  facts 
that  the  number  of  tests  required  to  obtain  the  necessary  experimental 
parameters  is  minimal,  and  the  generality  present  in  the  constitutive 
equations  is  adequate  to  describe  nonlinear  viscoelastic  behavior. 

In  summary,  while  the  empirical  approaches  are  useful,  they  do  not 
readily  allow  for  the  incorporation  of  temperature  or  moisture  effects. 

On  the  other  hand,  mechanical  models  are  very  difficult  if  not  impossible 
to  visualize  for  a  combined  stress  state.  Differential  operator  methods 
rely  upon  linear  solution  techniques,  which  restrict  the  applicable 
range  of  solutions.  Multiple  integral  representations  can  be  extremely 
difficult  to  evaluate  experimentally. 

As  will  be  shown  in  Chapter  IV,  it  is  quite  easy  to  incorporate 
the  measured  material  parameters  into  the  equations  derived  by  Schapery 
[8].  Also,  it  has  been  demonstrated  by  Schapery  [11]  that  for  a  multi- 
axially  stressed  material,  the  properties  can  be  expressed  as  a  function 
of  the  local  octahedral  shear  stress  in  the  material.  This  allows 
one  to  completely  characterize  a  material  by  a  series  of  uniaxial  tests. 
In  order  to  include  the  effects  of  temperature  and  moisture,  one  need 
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simply  repeat  the  above  uniaxial  tests  for  a  range  of  constant  temper¬ 
ature  and  moisture  levels,  to  determine  the  dependence  of  certain  con¬ 
stants  in  the  constitutive  equations  upon  these  parameters. 

A  detailed  description  of  the  analysis  method,  utilizing  the  con¬ 
stitutive  equations  developed  by  Schapery  [8]  in  a  finite  element  micro¬ 
mechanical  analysis,  is  given  in  the  following  chapter.  Succeeding 
chapters  present  some  numerical  results  and  conclusions  based  upon  these 


results. 


CHAPTER  III 


ANALYSIS  METHOD 


A  Statement  of  Linearity 

It  is  important  to  first  examine  the  concept  of  linear  viscoelastic 
material  behavior  and  its  range  of  applicability,  in  order  for  nonlinear 
viscoelastic  constitutive  equations  to  be  meaningful. 

While  certain  materials  may  behave  in  a  linear  manner  up  to  some 
relatively  high  stress  level  (for  a  fixed  temperature,  and  for  certain 
composites,  moisture  content),  others  may  be  quite  nonlinear  in  that 
same  stress  range.  In  order  to  differentiate  between  linear  and  non¬ 
linear  behavior,  a  definition  of  linearity  is  necessary. 

There  are  two  criteria  which  must  be  met  in  order  for  material  re¬ 
sponse  to  be  considered  linear  [1,12].  The  first  of  these  is  the  prop¬ 
erty  of  homogeneity.  This  property,  sometimes  referred  to  as  the  pro¬ 
portionality  property,  states  that  the  material  response  is  linearly 
proportional  to  the  input  history.  That  is,  if  the  stress  (input)  is 
doubled,  the  strain  response  (output)  is  correspondingly  doubled.  This 
can  be  written  as 

e{Ca(t) }  =  Ce{o(t) )  (2) 

where  C  is  any  constant,  o(t)  is  the  stress  input  function,  and  e{)  is 
the  material  response  based  upon  the  entire  history  of  o(t). 

The  second  criterion  is  the  superposition  principle.  This  property 
states  that  the  response  due  to  a  multistep  input  of  stress  with  time 
is  identical  to  the  sum  of  the  responses  due  to  individual,  independent 


stress  inputs  at  corresponding  tines.  This  can  ,be  written  as 

i (t)  +  02(t)  +  •••}  =  e(oi(t)}  +  e{o2(t)}  +  •••  (3) 

where  ojCt),  a2(t),  •••  represent  any  type  of  stress  input  history. 

This  second  property  can  be  demonstrated  by  an  example.  Assume 
three  identical  pieces  of  material,  each  subjected  to  a  different  load¬ 
ing  history.  Specimen  A  is  subjected  to  a  load  input  described  by 
Figure  2a.  Specimen  B  by  that  of  Figure  2b,  and  Specimen  C  by  that  of 
Figure  2c.  If  the  material  behavior  is  linear  viscoelastic,  the  super¬ 
position  principle  states  that  the  sum  of  the  strain  responses  from  the 
loadings  shown  in  Figures  2b  and  2c  will  be  identical  to  the  response 
from  the  loading  shown  in  Figure  2a. 
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a)  Stress  input 
history  for 
Specimen  A 

Figure  2.  Plots  of 


b)  Stress  input 
history  for 
Specimen  B 

versus  time  for  three 


c)  Stress  input 
history  for 
Specimen  C 

different  specimens 


Many  investigators  assume  linearity  if  a  material  satisfies  the 
hotnogenlty  principle  alone,  while  in  fact  the  material  may  be  quite  non¬ 
linear.  It  can  be  shown  [1,12]  that  if  the  material  satisfies  the  super- 


position  criterion,  it  will  automatically  satisfy  homogenity,  but  the 
opposite  is  not  true.  Thus,  if  a  material  is  to  be  examined  for  linear 
viscoelastic  behavior,  the  superposition  test  must  be  applied. 


Constitutive  Equations 

The  single  integral  nonlinear  viscoelastic  constitutive  equation 
derived  by  Schapery  [8]  for  an  isotropic  material  under  isothermal  con¬ 
ditions  and  uniaxial  loading  can  be  written  as 
e(t)  =  g0D0a(t)  + 

g!  (  ADOF-V’)  dtr  (4) 

J  0~  dr' 

where  D0  is  the  initial  (time-independent)  value  of  the  compliance  and 
ADCV)  is  the  transient  (time-dependent)  value  of  the  compliance.  Note 
that  AD(T)  is  a  function  of  reduced  time  Y.  The  reduced  times  are  defined 
as , 

V  =  f  —  and  r  = 

Jo  aq 

The  terms  go,  gj ,  g2,  and  Sq  are  all  functions  of  stress  and  have  cer¬ 
tain  thermodynamic  significance  [8];  for  example,  gg,  gi ,  and  g2  reflect 
stress-dependence  on  Gibb's  free  energy  while  changes  in  a0  signify 
effects  of  both  entropy  production  and  free  energy.  The  term  a0  in  Eq. 

(5)  is  also  a  function  of  temperature,  and  for  certain  composites,  mois¬ 
ture  content.  This  factor  tends  to  accelerate  or  decelerate  the  influence 
of  time,  meaning  that  changes  in  temperature,  moisture,  or  stress  produce 
the  same  output  of  events  (with  respect  to  some  reference  value),  but  at 
a  faster  or  slower  rate. 
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Multiaxial  relationships  as  defined  by  Schapery  [8]  are  somewhat 
more  complicated  in  that  they  allow  for  any  degree  of  coupling  between 
the  applied  stresses.  In  the  present  analysis  it  is  assumed  that  this 
coupling  (except  Poisson  effects)  is  of  higher  order  and  is  neglected. 

The  constitutive  equations  used  in  the  present  analysis  can  be  written 
in  the  simple  form 

£ij(t>  -  Di;]kl(t ,T,M, a)okl  (6) 

where  Djjk^(t,T,M,o)  is  the  time-dependent  nonlinear  creep  compliance 
tensor,  and  is  a  function  of  time,  temperature,  moisture,  and  applied 
stress  level. 

In  order  to  useEq.  (6)  in  a  finite  element  analysis  scheme,  it  be¬ 
comes  necessary  to  simplify  the  convolution  integrals  (Eq.  4)  contained 
in  D±jkl. 

In  the  current  analysis  an  incremental  procedure  is  used,  all  load¬ 
ing,  temperature,  and  moisture  increments  being  applied  in  time-independ¬ 
ent  steps.  That  is,  these  increments  are  assumed  to  occur  instantaneous¬ 
ly,  and  to  be  held  constant  during  succeeding  time  steps.  Time-depend¬ 
ent  response  is  determined  by  the  application  of  finite  increments  of 
time;  no  increments  of  stress,  temperature,  or  moisture  are  permitted 
during  a  time  increment.  In  other  words,  stress  is  a  constant  during 
any  time  increment. 

By  substituting  a(t)  =  aa  (a  constant  value)  and  integrating,  Eq. 

(4)  yields 

e(t)  -  [g0D0  +  g1g24D(t/a0)]oa.  (7) 

Details  of  this  integration  are  presented  in  Appendix  A.  Note  that  AD 
is  a  function  of  the  ratio  of  real  time  to  the  time  shift  factor  ag. 
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From  Eq.  (7),  the  nonlinear  creep  compliance  for  a  constant  applied 
stress  oa  is  defined  as 

D  -  =  g0Do  +  8182  AD(t/a0)  (8) 

where  the  term  go^o  i8  the  instantaneous  (time- independent)  portion  and 
the  term  gig2  AD(t/aa)  is  the  transient  (time-dependent)  portion  of  the 
creep  compliance. 

As  mentioned  previously,  no  coupling  between  multiaxial  loading  is 
assumed  except  for  Poisson  effects.  If  the  material  is  isotropic  (as  is 
typically  the  case  of  the  matrix  material  in  a  composite) ,  the  creep 
compliance  tensor  in  Eq.  (6)  is  simplified  to  three  time-dependent  terms, 
two  of  which  are  independent. 

The  isotropic  elastic  constitutive  equations  for  multiaxial  loading, 
written  in  tensorial  form,  are  [13] 

£ij  =  °ij  "  I  akk6iJ  (9) 

where  E  is  the  elastic  (Young's)  modulus  and  v  is  Poisson's  ratio.  By 
using  Eq.  (8)  (for  the  special  case  of  constant  externally  applied 
stresses),  making  Poisson's  ratio  time-dependent,  and  substituting  into 
Eq.  (9),  the  time-dependent  constitutive  equations  (Eqs.  6)  become 

tjj (t)  =  D[1  +  v(t)]oij  -  Dv(t)akkeij  (10) 

Eqs.  (10)  are  valid  for  any  isotropic,  time-dependent  material  subject 
to  a  constant  externally  applied  load. 

Method  of  Solution 

In  the  present  study,  it  was  desired  to  model  a  unidirectional  com¬ 
posite  material  by  using  a  finite  element  micromechanics  analysis.  A 


finite  element  grid  of  the  internal  structure  of  a  composite  is  construc¬ 
ted,  allowing  the  examination  of  the  internal  stresses  in  a  composite 
due  to  a  variety  of  external  conditions.  The  cross-sectional  arrange¬ 
ment  of  fibers  in  most  unidirectional  fibrous  composites  typically  tends 
to  be  random.  However,  there  is  usually  little  loss  of  generality  in 
the  analysis  if  some  periodicity  is  assumed,  as  shown  in  Figure  3. 


Figure  3.  Fiber  packing  arrangement  of  a  unidirectional  composite 

A  detailed  discussion  of  the  various  possible  types  of  fiber  packing,  and 
the  ramifications  of  assuming  each  type,  are  presented  by  Miller  and  Adams 
[14]  and  will  not  be  repeated  here.  In  the  present  analysis,  the  fibers 
are  assumed  to  be  arranged  in  a  rectangular  array,  as  shown  in  Figure  4a. 
Since  this  particular  type  of  packing  is  assumed,  symmetry  arguments  can 

be  used  [14],  requiring  only  one  quarter  of  a  fiber  to  be  modeled,  as 

\ 

shown  in  Figure  4b.  \ 

In  order  to  simulate  a  small  continuum  in  a  rectangular  array  of 
fibers  in  a  composite,  certain  boundary  conditions  must  be  imposed. 


* 


Figure  4.  Unit  cell  for  the  finite  element  analysis 


For  descriptive  purposes,  a  coordinate  system  has  been  set  up  in  Figure 
4b  with  the  origin  located  in  the  lower  left  hand  corner.  The  right  hand 
boundary  corresponding  to  x»a  is  assumed  to  displace  parallel  to  the 
original  boundary.  That  is,  all  displacements  in  the  x-direction  along 
the  boundary  x=a  are  equal  in  magnitude.  Similarly,  all  displacements 
in  the  y-direction  along  the  boundary  y=b  are  considered  to  be  equal. 

Both  of  these  boundaries,  x=a  and  y=b,  are  allowed  to  move  freely  in  the 
y-  and  x-directions,  respectively.  The  remaining  two  boundaries,  x=0 
and  y=0,  are  assumed  to  be  fixed  in  the  x-  and  y-directions ,  respectively; 
yet  they  are  free  to  move  in  the  y-  and  x-directions,  respectively.  Under 
the  conditions  of  generalized  plane  strain,  which  will  be  discussed 
shortly,  all  out-of-plane  normal  strains  ez  are  assumed  to  be  the  same. 

In  other  words,  all  displacements  in  the  z-direction  are  equal. 

The  finite  element  grid  used  in  the  present  analysis  is  shown  in 
Figure  5.  A  detailed  description  and  flow  chart  of  the  computer  program 
(including  time-dependence)  is  presented  in  Appendix  B. 

In  the  present  micromechanics  analysis,  thermal-  and  moisture- 
dilatational  effects  are  included  by  means  of  initial  strains  [15].  In 
order  for  the  creep  formulation  to  be  compatable  with  the  current  pro¬ 
gram,  a  scheme  including  creep  strains  as  initital  strains  was  adopted. 

The  technique  for  accomplishing  this  will  be  briefly  explained  in  the 
following  paragraphs  and  illustrated  in  detail  by  mathematical  equations 
after  that. 

First,  however,  it  is  appropriate  to  convert  Eq.  (10),  which  is  in 
tensorial  form,  to  the  notation  used  by  Miller  and  Adams  [14].  The 
time-dependent  constitutive  equation  for  an  isotropic  material  under  a 
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The  finite  element  grid  used  in  the  present  analysis 
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generalized  plane  strain  formulation  for  a  constant  stress  state  is 


r  n 

(  > 

£x(t) 

1  — v (t )  0  — v(t) 

<*x 

£y(t) 

1  0  -v(t) 

CTy 

)  =  D 

l 

YxyCt) 

l+v(t)  0 

1 

Txy 

ez(t) 

symmetric  1 

°z 

V  J 

Under  conditions  of  generalized  plane  strain,  the  normal  strain  e2  is  a 
constant  over  the  entire  model  (Figure  5),  as  previously  mentioned,  and 
the  shear  strains  yx2  and  Yyz  are  assumed  zero.  In  Eq.  (11),  D  is  given 
by  Eq.  (8)  and  v(t)  is  the  time-dependent  Poisson's  ratio. 

In  the  finite  element  analysis,  nodal  point  forces  are  generated 
due  to  the  creep  strains.  The  magnitudes  of  these  induced  nodal  forces 
are  a  function  of  the  geometry  of  an  element  and  the  time-dependent  ma¬ 
terial  properties  of  that  element.  These  node  point  forces  are  calculat¬ 
ed  using  Eq.  (11)  and  the  elemental  stresses  at  the  beginning  of  the 
time  increment.  Since  this  is  an  incremental  finite  element  program, 
care  must  be  taken  to  insure  that  only  the  incremental  and  not  the  total 
creep  strains  are  calculated.  To  accomplish  this,  Eq.  (11)  is  differen¬ 
tiated  with  respect  to  time  and  evaluated  at  the  time  corresponding  to 
the  midpoint  of  the  time  increment.  By  assuming  the  elemental  stresses 
to  remain  constant  throughout  the  time  interval,  and  multiplying  through 
by  the  time  increment  At,  a  linear  approximation  of  the  incremental  creep 
strain  for  that  time  increment  is  obtained.  For  this  reason  it  is  nec¬ 
essary  to  keep  the  time  increments,  which  are  input  into  the  analysis, 
very  small  when  the  value  of  the  creep  compliance  is  changing  rapidly. 

For  a  constant  externally  applied  stress,  this  rapid  change  corresponds 
to  the  initial  (or  primary)  part  of  the  creep  curve.  Thus,  it  is 
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necessary  to  keep  the  time  increments  small  just  after  a  load  has  been 
applied  or  incremented.  The  magnitudes  of  these  initial  time  increments 
will  of  course  depend  upon  the  viscoelastic  properties  or  the  type  of 
material  being  modeled.  Polymers  may  require  many  small  increments  (of 
minutes  or  fractions  of  minutes) ,  while  steels  may  require  only  a  few 
increments  (of  minutes) . 

Details  of  the  strain-displacement  relations,  shape  functions,  and 
basic  finite  element  formulation  for  the  time-independent  formulation  are 
given  by  Miller  and  Adams  [14]  and  will  not  be  repeated  here.  Rather, 
only  the  equations  which  deal  with  the  time-dependent  finite  element 
formulation  will  be  presented. 

The  equation  which  relates  nodal  forces  to  displacements  can  be 
written  as  [14] 

•  •  • 

{F}i  =  [K^CSJi  +  {Fe  }±  (12) 

where  {F}^  is  the  incremental  nodal  force  vector  of  the  it*1  element  due 

to  externally  applied  loads,  (F_  is  the  incremental  nod*l  force  vec- 

e0 

tor  due  to  temperature,  moisture,  and  time  changes,  [K]i  is  the  element 
stiffness  matrix,  and  {6}^  is  the  incremental  element  nodal  displacement 
vector. 


During  a  time  increment,  {F}^  must  be  zero  since  loads  can  only  be 
applied  during  time- independent  increments.  Also,  no  changes  in  temper¬ 


ature  or  moisture  content  can  occur.  Thus,  {Fe  is  the  incremental 


nodal  force  due  to  incremental  strain  changes  caused  by  creep  or  the 


time-dependent  properties  only. 

The  element  stiffness  matrix  is  defined  as  [14] 
[K] J  -  [B]J[D]1[B]it1A1 


(13) 


where  [B]^,  sometimes  referred  to  as  the  shape  matrix,  relates  the  node 
point  displacements  and  the  element  strains,  t^  is  the  thickness  of  the 
element  in  the  out-of-plane  direction  (which  is  unity  in  the  present  an¬ 
alysis),  and  is  the  cross-sectional  area  of  the  element.  The  element¬ 
al  material  properties  matrix  [D]j[  is  found  by  evaluating  the  time-depend¬ 
ent  parameters  such  as  D  (see  Eq.  8)  and  v(t)  at  the  time  corresponding 
to  the  midpoint  of  the  time  increment,  and  then  substituting  into  the 
following  matrix  for  an  isotropic  material  under  generalized  plane  strain 


conditions : 


1  v*/(l-v*)  0  v*/(l-v*) 

1  0  v*/(l-v*) 


[D]  = 


(1-v*) 


D*(l+\j*)  (l-2v*) 


[(l-2v*)/2(l-v*) ]  0 


symmetric 


The  terms  D*  and  v*  are  the  creep  compliance  and  the  time-dependent 
Poisson's  ratio,  respectively,  evaluated  at  the  particular  time  defined 


above. 


Using  Eq.  (12)  and  the  finite  element  analysis,  the  incremental 

displacements  can  now  be  solved  for  if  the  incremental  forces  {F^}^  due 

to  the  creep  strains  are  known.  The  expression  for  (Fe  L  is  [14] 

0  i 

(F^i  =  [B]J  [D]{eo>t1A1  (15) 

where  {eQ}  is  the  incremental  creep  strain. 

Evaluation  of  {eQ}  for  an  element  is  accomplished  by  the  use  of 
Eqs.  (11)  and  the  element  stresses  at  the  beginning  of  the  increment. 

To  get  {eQ},  the  derivative  with  respect  to  time  is  taken  of  Eqs.  (11) 
and  all  of  the  time-dependent  properties  evaluated  at  the  time  corres¬ 
ponding  to  the  midpoint  of  the  time  interval.  Using  the  first  of  Eqs. 
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(11)  as  an  example,  i.e., 

ex  =  Dox  -  Dv(t)Oy  -  Dv(t)oz 
its  derivative  with  respect  to  time  is 

|j(ex)  -  -  (t^(n)lv(c)  + 

-  ([^(D)Mt)  +  [|jv<t))D}oz 


(16) 


(17) 


It  has  been  determined  by  several  investigators  [12,16,17]  that  Poisson's 
ratio  is  time-independent  or  a  very  weak  function  of  time.  For  this 
reason,  it  was  decided  that  certain  terms  in  Eq.  (17)  could  be  considered 
higher  order  effects  and  therefore  be  dropped  from  the  equations,  viz.,  , 
the  terms  involving  time  rates  of  change  of  Poisson's  ratio.  With  this 
simplification,  Eq.  (17)  becomes 


ft<ex>  “  ^(°)ox  -  &(D)]v(t)oy  -  [fF(D)]v(t)oz.  (18) 

The  terms  gQ,  DQ,  gj ,  g2,  and  a0  in  the  nonlinear  creep  compliance  D 
(Eq.  8),  are  functions  of  stress,  temperature  and  moisture  only,  and 
are  required  to  be  constant  throughout  a  time  increment  in  this  analysis. 
Thus,  when  Eq.  (18)  is  evaluated  at  a  time  corresponding  to  the  midpoint 
of  the  time  increment  and  multiplied  through  by  At,  a  linear  approxima¬ 
tion  of  the  incremental  creep  strain  ex  is  obtained. 

A  similar  argument  is  applied  to  the  relations  for  cy,  ez,  and  YXy; 
the  resulting  equations  are  given  below. 


JJ-(sy)  -  -  l^(i)Jv(t)ox  +^(5)oy  -  [^<WMc)oz 

-  -  i|f(D)]v(t)ox  -  [^(D)]v<t)oy  +  5j(D)„2  (19) 

Jj(yxy>  -  &(i>][i+v<t>]Txy 


Once  these  incremental  values  of  strain  are  known,  the  program 
solves  for  the  incremental  node  point  displacements,  which  are  later 
used  to  find  the  new  values  for  elemental  stress.  When  this  has  been 
accomplished,  the  program  returns  for  the  next  increment  of  time,  tem¬ 
perature,  moisture,  or  load. 

Multiple  Step  Loading 

For  a  two-step  loading,  as  shown  in  Figure  6,  the 


Figure  6.  Plot  of  stress  versus  time  for  a  two-step  input 


constitutive  equation  is  given  by 


e(t)  =  gQD0ob  +  g^o[ 


Pb 


b  a._ 
g  g  AD 
1  2 


t-ta 


°a  + 


gbgaADj-I£.  + 
1  2 


t-t. 


(20) 


for  time  t>ta.  Details  of  the  derivation  of  Eq.  (20)  are  presented  by 
Schapery  [18]  and  also  in  Appendix  A.  The  superscripts  and  subscripts 
a  and  b  are  used  to  denote  which  constants,  stresses,  and  times  are 
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associated  with  each  stress  level.  It  is  important  to  note  the  two 
additional  terms  which  appear  in  Eq.  (20),  that  did  not  appear  in  Eq.  (7). 
These  two  terms  represent  the  application  of  a  -c  stress  and  a  +o^ 
stress,  beginning  at  time  t=ta.  It  is  necessary  to  include  the  two  ad¬ 
ditional  terms  due  to  the  nonlinear  dependence  of  certain  constants  in 
Eq.  (8)  upon  the  applied  stress  level.  If  the  material  behavior  were 
linearly  viscoelastic,  it  would  be  necessary  to  add  only  one  term  to 
Eq.  (7),  viz.,  a  term  involving  a  stress  increment  of  the  difference 
between  the  two  stress  inputs  (ob-0a^ • 

While  other  authors  [19,20]  have  developed  recursive  relationships 
which  require  integration  of  the  viscoelastic  constitutive  relations  over 
the  current  time  increment  only,  no  feasable  method  of  keeping  track  of 
the  stress  history  was  discovered  which  was  compatable  with  the  current 
finite  element  formulation.  This  means  that  individual  histories  for 
each  load  step  had  to  be  kept  track  of  in  the  current  analysis  and  the 
contribution  of  each  history  solved  for  separately.  For  example,  a  two- 
step  loading  shown  in  Figure  6  requires  the  program  to  solve  for  three 
contributions  of  the  total  creep  strain  for  time  t>ta;  the  strain  due  to 
the  initial  application  of  +ca  at  time  t=0,  the  strain  due  to  a  pseudo 
application  of  -oa  at  time  t=ta,  and  the  strain  due  to  a  pseudo  applica¬ 
tion  of  +Oj,  at  time  t=ta.  Each  additional  step  in  stress  would  require 
two  extra  solutions.  While  the  program  is  capable  of  any  number  of  step 
loadings,  it  can  realistically  model  only  a  finite  number  of  step  loadings 
due  to  the  increase  of  computer  time  necessary  to  solve  additional  in¬ 
crements.  Using  a  CDC  Cyber  760  computer,  approximately  2-3  seconds  of 
computer  time  are  required  to  solve  one  stress  increment  using  the 
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finite  element  mesh  shown  in  Figure  5.  (This  mesh  has  a  bandwidth  of 
75.)  Thus,  for  times  t>ta,  approximately  6-9  seconds  of  computer  time 
are  required  for  each  time  increment. 

While  the  formulation  by  Haisler  and  Saunders  [19,20]  avoids  this 
difficulty,  it  requires  the  use  of  an  iterative  scheme  and  hence  cannot 
readily  model  materials  with  "flat"  stress-strain  curves  such  as  aluminum 
[21].  Thus,  the  micromechanical  modeling  of  a  composite  such  as  boron/ 
aluminum  would  be  very  difficult.  Also,  iterative  schemes  can  require 
large  amounts  of  computer  time  to  converge  upon  a  solution. 

The  intent  of  the  present  study  was  not  to  analyze  the  relative  ad¬ 
vantages  and  disadvantages  of  a  particular  finite  element  method.  Rather, 
the  goal  was  to  use  an  existing  finite  element  micromechanics  program  and 
incorporate  time-dependent  effects.  In  this  chapter,  the  viscoelastic 
constitutive  equations  were  presented  and  a  scheme  for  incorporating 
these  equations  into  a  finite  element  analysis  was  devised.  However, 
numerical  values  for  the  time-dependent  parameters  must  be  known.  In 
the  following  chapter,  a  brief  experimental  program  is  described  which 
was  conducted  to  determine  these  parameters. 


CHAPTER  IV 


EXPERIMENTAL  CHARACTERIZATION 

In  order  to  utilize  the  constitutive  equations  presented  in  the 
previous  chapter,  a  series  of  creep-recovery  tests  is  necessary.  To  in¬ 
clude  the  dependence  of  the  kernels  in  the  equations  upon  temperature 
and  moisture,  a  full  experimental  program  of  isothermal,  constant  humid¬ 
ity,  and  constant  stress  tests  must  be  conducted.  A  method  of  evaluating 
all  parameters  given  in  Eqs.  (4)  and  (5)  has  been  developed  by  Lou  and 
Schapery  [11]  and  will  be  demonstrated  in  this  chapter. 

Test  Method 

In  the  present  study,  a  brief  experimental  program  was  undertaken 
in  order  to  evaluate  the  creep  constants  which  were  needed  as  input  to 
the  finite  element  program.  This  experimental  program  was  not  intended 
to  completely  evaluate  the  viscoelastic  parameters  necessary  for  the  con¬ 
stitutive  equations.  Rather,  the  goal  was  to  obtain  approximations  of 
the  numerical  values  for  the  kernels  to  be  used  as  input  to  the  computer 
program.  A  series  of  compression  creep-recovery  tests  were  run  on 
Hercules  3501-6  epoxy  resin  [22]  test  specimens.  A  compression  test  was 
selected  to  avoid  the  gripping  problem  usually  associated  with  tensile 
tests  of  neat  (unreinforced)  resin  materials  [121. 

For  this  program,  a  total  of  fifteen  specimens  were  tested  at  five 
different  stress  levels  (i.e.,  three  specimens  at  each  stress  level 
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selected).  All  tests  were  conducted  at  room  temperature  (21°C)  and 
ambient  humidity. 

The  specimens  were  fabricated  from  *2- inch  diameter  rods  cast  in  a 
split  steel  mold.  The  rods,  molded  in  approximately  6-inch  lengths, 
were  cut  into  1*2- inch  long  test  specimens.  The  ends  were  then  ground 
flat  using  a  surface  grinder  to  insure  that  they  were  smooth  and  parallel 
to  each  other.  Details  of  the  fabrication  and  curing  process  may  be 
found  in  Reference  [23]. 

An  Instron  Model  1321  servohydraulic  testing  machine  was  used  to 
apply  and  maintain  a  constant  load  on  the  rods.  Strain  was  measured 
using  a  %-inch  gage  length  extensometer  and  all  data  were  taken  using  a 
Hewlett-Packard  21MX  Series  E  minicomputer.  The  values  of  load,  strain, 
and  time  were  taken  just  prior  to  and  just  after  application  of  the  load. 
The  total  time  to  apply  load  for  all  tests  was  about  four  seconds.  After 
application  of  the  load,  values  of  load,  strain,  and  time  were  recorded 
each  minute  for  ten  minutes.  After  that  time,  data  were  recorded  every 
five  minutes.  The  creep  tests  lasted  for  approximately  one  hour,  after 
which  time  the  computer  automatically  removed  the  load  while  continuing 
to  record  the  load,  strain,  and  time.  At  this  point  data  were  taken  at 
one  minute  intervals  for  ten  minutes.  After  ten  minutes,  data  were 
taken  every  five  minutes  for  another  two  hours. 

This  procedure  was  used  for  each  of  five  different  compressive 
stress  levels,  i.e.,  2000,  6000,  8000,  16,000,  and  20,000  psi.  The 
tests  conducted  at  2000  psi  were  checked  for  linearity  and  found  to 
satisfy  the  principle  of  superposition  as  discussed  in  Chapter  III. 
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Experimental  Results 

Plots  of  the  data  taken  are  shown  in  Appendix  C.  These  plots  are 
actual  unretouched  data  acquired  by  the  minicomputer. 

The  technique  developed  by  Lou  and  Schapery  [11] ,  which  reduces  the 
data  presented  in  Appendix  C,  makes  use  of  a  power  law  representation  to 
describe  the  nonlinear  creep  compliance  expressed  by  Eq.  (8).  This  meth¬ 
od  utilizes  recovery  data  to  determine  the  exponent  used  in  the  power 
law.  A  constitutive  equation  describing  the  recovery  portion  of  the 
creep  curve  is  derived  in  Appendix  A  and  the  equations  are  repeated  here 
for  convenience.  For.  the  stress  history  described  in  Figure  7,  the  re¬ 
covery  strain  is 


where 


er(t) 


^  [(1  +  aaA)n  -  (aaA)n] 


(21) 

(22) 


and 

Aea  =  e(ta)  -  eQ  =  g^Ct^  (23) 

for  time  t>ta.  The  term  C  is  a  constant  from  the  power  law  (see  Appendix 
A)  and  the  term  eQ  in  Eq.  (23)  is  the  initial  or  time-independent  value 
of  strain  which  occurs  at  time  t=0. 


Figure  7.  Plot  of  stress  versus  time  for  a  creep-recovery  test 
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The  first  step  in  the  evaluation  of  the  creep  parameters  is  to  de¬ 
termine  the  exponent  n.  This  may  be  found  by  plotting  the  equation 

er(t)  =  (1  +X)n  -  (X)n  (24) 

on  log-log  paper  for  several  values  of  n  (usually  0<n<0.5)  and  overlay¬ 
ing  these  curves  on  the  experimental  recovery  curves.  Equation  (24)  is 

generated  by  setting  Aea,  g  ,  and  equal  to  one  in  Eq.  (21).  This  is 

1  0 

done  since  only  the  exponent  is  being  determined  at  this  time.  The  best 
fit  to  all  of  the  data  will  determine  a  value  for  n. 

Since  the  tests  run  at  2000  psi  satisfied  the  superposition  principle, 
the  response  was  linear  viscoelastic  and  all  of  the  terms  gQ,  g^,  and 
a0  were  equal  to  one.  For  this  reascn,  these  particular  curves  were  aver¬ 
aged  and  used  as  a  reference  curve  for  all  other  stress  levels  tested. 

By  using  a  reference  curve,  the  values  of  g^,  g^ ,  g^,  and  a„  as  a  function 
of  stress  were  determined. 

) 

To  make  the  recovery  curves  described  by  Eq.  (21)  match  the  refer¬ 
ence  curve  (at  2000  psi) ,  a  certain  amount  of  horizontal  and  vertical 
shifting  is  necessary.  The  amount  of  horizontal  shifting  required  deter¬ 
mines  the  value  for  a  .  If  a0>l  the  data  must  be  shifted  to  the  right, 
and  if  aQ<l  the  data  must  be  shifted  to  the  left.  The  amount  of  vertical 
shifting  determines  the  value  of  Aca/gj .  If  Aea/gj>l  the  data  must  be 
shifted  downward;  if  Aea/gj<l  the  data  must  be  shifted  upward. 

Once  the  values  for  a0  and  Aea/gj  have  been  determined,  the  value  of 
gj  may  be  calculated  using  the  first  part  of  Eq.  (23)  to  get  Aea.  Plots 
of  a0  and  gj  versus  stress  for  Hercules  3501-6  epoxy  resin  are  shown  in 
Figures  8  and  9,  respectively.  Note  that  these  plots  show  abrupt  changes 
around  2000  psi.  This  was  intentionally  done  to  force  these  constants  to 


versus  stress  for  Hercules  3501-6  epoxy  resin 


have  a  value  of  one,  during  the  range  of  linear  viscoelastic  behavior. 

It  would  not  be  meaningful  for  these  constants  to  have  values  other  than 
one  in  this  range.  To  determine  exactly  how  the  curves  transition  from 
linear  to  nonlinear  viscoelasticity,  more  tests  need  to  be  run  at  these 
stress  levels. 

The  only  remaining  creep  constants  which  need  to  be  determined  are 
gQ  and  g2.  Values  of  the  product  g0DQin  Eq.  (8)  are  easily  calculated 
since  they  are  merely  the  ratio  of  initial  strain  to  the  applied  stress 
(eQ/oa).  A  plot  of  the  relationship  of  the  product  gQD0  to  stress  is 
shown  in  Figure  10. 

The  functional  relationship  of  g 2  with  stress  was  determined  by 
utilizing  a  power  law  in  Eq.  (7).  Substitution  of  the  power  law  Eq. 

(A-15)  into  Eq.  (7)  gives 

e(t)  =  g0D0aa  +  cg1g2(^-]  oa.  (25) 

\  °l 

Using  Eq.  (25)  and  the  parameters  ggDg,  gj ,  and  a0  just  determined  from 
the  creep  data,  the  product  Cg2  may  be  determined  for  each  value  of  stress. 
During  the  linear  viscoelastic  range  (up  to  approximately  2000  psi  in  the 
present  case) ,  g2  is  equal  to  one  and  C  is  equal  to  the  net  time-depend¬ 
ent  (transient)  portion  of  the  creep  compliance  at  unit  time.  By  divid¬ 
ing  through  by  C  for  each  product  Cg2  determined,  the  term  g2  as  a 
function  of  stress  is  obtained.  A  plot  of  g2  versus  stress  for  Hercules 
3501-6  epoxy  resin  is  as  shown  in  Figure  11. 

In  Appendix  C,  Figures  C-ll  through  C-15  show  actual  data  (indicated 
by  the  stars)  plotted  against  the  predicted  creep  response,  using  the 
creep  parameters  just  calculated.  There  is  excellent  agreement  between 
the  data  and  the  prediction,  as  expected. 
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Some  difficulty  was  encountered  when  shifting  of  the  recovery  data 
was  done  to  obtain  the  terms  aa  and  Aea/gj.  For  this  reason,  the  data 
plotted  in  Figures  8  through  11  may  be  questionable.  In  particular,  the 
amount  of  shifting  done  for  each  test  was  not  always  consistent.  For 
some  of  the  tests,  a  specific  amount  of  horizontal  and  vertical  shifting 
was  done,  while  for  others  (subjected  to  the  same  stress  level)  only  ver¬ 
tical  shifting  was  necessary.  Thus,  the  results  plotted  in  Figures  8  and 
9  have  been  obtained  through  a  degree  of  estimation.  Agreement  between 
the  data  and  the  predictions  presented  in  Appendix  C  for  creep  is  good 
since  the  final  parameter  g 2  is  found  by  using  creep  data  and  Eq.  (25), 
with  the  values  determined  for  DQg0,  gj ,  and  a0.  The  problem  of  a  poor 
fit  of  the  predictions  to  the  data  would  arise  only  when  the  kernels 
appear  individually  and  not  together  as  the  product  g1g2/(a0)n  in  Eq.  (25). 
One  such  instance  would  be  during  recovery.  In  Eq.  (21),  the  kernel  gj 
cancels  with  gj  contained  in  the  term  Aea  and  hence  does  not  appear  in 
the  equation.  This  problem  will  be  exemplified  in  the  next  chapter  by 
recovery  predictions  made  using  these  kernels.  To  avoid  this  difficulty, 
data  should  be  sampled  more  often  during  both  the  creep  and  recovery  tests. 
Also,  it  has  been  recommended  in  the  literature  [11,12]  that  a  curve  be 
drawn  through  the  data  points  and  used  for  shifting  instead  of  using  the 
actual  data  points  as  was  done  here.  This  helps  to  smooth  out  and  also 
eliminate  extaneous  data  points.  By  following  these  suggestions  it  is 
expected  that  more  accurate  values  for  gQ,  gj,  g2,  and  a0  will  be  obtained. 


CHAPTER  V 


NUMERICAL  EXAMPLES 

Several  examples  are  presented  in  this  chapter  to  demonstrate  the 
capability  of  the  analysis.  A  short  description  of  each  load  case  is  pre¬ 
sented,  followed  by  a  series  of  contour  plots  illustrating  the  internal 
stress  state  of  the  composite.  All  examples  use  the  finite  element  mesh 
shown  in  Figure  5,  which  corresponds  to  a  63  percent  fiber  volume,  except 
where  noted. 

The  computer  program  was  first  checked  out  by  applying  a  uniaxial 
tensile  load  to  a  model  made  up  entirely  of  resin.  Two  independent  sets 
of  data,  one  generated  in  this  study  for  Hercules  3501-6  epoxy  resin  and 
another  generated  by  Beckwith  [12]  for  Shell  58-68R  epoxy  resin,  were 
used  as  input  to  the  program.  Temperature  and  moisture  levels  were  held 
constant  during  each  simulation  so  that  the  solution  generated  by  the 
computer  program  could  be  compared  to  the  input  data. 

The  next  check  was  performed  by  modeling  an  actual  composite  sub¬ 
jected  to  a  uniaxial  transverse  normal  load.  The  results  predicted  by 
the  program  were  compared  with  actual  data  generated  by  Beckwith  [12] 
and  Irion  [24]  using  the  Shell  and  Hercules  resin  systems,  respectively. 

Once  verifications  with  simple  uniaxial  creep  tests  were  demonstrat¬ 
ed,  the  effect  of  fiber  volume  was  briefly  examined.  Next,  several  runs 
were  made  to  demonstrate  that  the  internal  stress  behavior  of  a  composite 
during  cooldown  can  be  predicted.  Finally,  uniaxial  cycling  and  multi- 
axial  load  tests  were  simulated  and  the  internal  stress  states  of  each 


composite  as  a  function  of  time  were  examined. 


The  analysis  is  capable  of  handling  the  effects  of  moisture;  however, 
no  examples  are  presented  in  this  chapter  since  there  were  no  experi¬ 
mental  data  available  concerning  this  effect.  The  absorption  of  moisture 
into  a  composite  has  the  same  general  effect  as  a  temperature  increase, 
i.e.,  both  have  a  similar  influence  on  the  time  shift  factor  aa. 

Specifically,  the  absorption  of  moisture  tends  to  make  a  composite  creep 
more.  As  mentioned  in  Chapter  IV,  to  include  the  effects  of  moisture, 
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several  creep-recovery  tests  at  various  moisture  levels  for  several  values 
of  stress  must  be  run. 

Creep  of  Epoxy  Resin 

Using  the  data  generated  by  Beckwith  [12]  for  the  Shell  epoxy  resin, 
a  uniaxial  load  was  applied  to  a  neat  epoxy  resin  model.  This  model  was 
made  by  specifying  all  elements  in  the  finite  element  mesh  shown  in 
Figure  5  to  have  neat  resin  properties.  A  stress  of  1000  psi  in  the  x- 
direction  (see  Figure  4b)  was  applied  for  2.1  hours  and  then  removed, 
after  which  recovery  was  determined  for  another  2.1  hours.  This  pro¬ 
cedure  was  repeated  for  three  separate  constant  temperature  levels,  i.e., 
-6.7°C  (20°F) ,  23.9°C  (75°F) ,  and  60.0°C  (140°F).  During  each  simulation, 
the  moisture  content  was  assumed  to  be  zero. 

Results  of  the  comparison  are  shown  in  Figures  12  and  13.  The  actual 
input  response  is  represented  by  the  solid  lines  while  the  predictions 
made  using  the  finite  element  analysis  are  shown  by  the  symbols.  To 
obtain  the  solid  lines  or  the  input  (expected)  response,  the  consti¬ 
tutive  equation  for  each  test  condition  was  written  and  then  plotted  for 


Figure  12.  Plots  of  strain  versus  time  for  Shell  58-68R  epoxy  resin  subject  to  a  2.1  hour 
creep  and  a  2.1  hour  recovery  test.  Predicted  values  are  shown  by  symbols 
while  the  solid  lines  represent  the  actual  experimental  response  [12]. 
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Figure  13.  Expanded  scale  plots  of  strain  versus  time  for  Shell  58-68R  epoxy  resin  subject  to  i 
2.1  hour  creep  and  a  2.1  hour  recovery  test.  Predicted  values  are  shown  by  symbols 
while  the  solid  lines  represent  the  actual  experimental  response  [12]. 


a  total  time  of  five  hours.  This  was  easily  done  since  the  moisture, 
temperature,  and  stress  levels  were  held  constant  throughout  the  tests 
(except  during  recovery,  where  the  stress  was  removed).  In  Figure  13 
the  strain  scale  has  been  expanded,  greatly  exaggerating  the  error  be¬ 
tween  the  predictions  (symbols)  and  the  expected  response  (solid  lines) . 
By  closely  examining  Figure  13,  it  can  be  seen  that  the  agreement  between 
the  predictions  and  the  expected  response  is  excellent;  the  maximum  error 
for  the  highest  temperature  (60.0°C)  is  less  than  3  percent.  To  obtain 
this  close  agreement  at  the  high  temperature,  several  time  increments  of 
one  minute  were  used  at  the  beginning  of  the  creep  and  recovery  portions 
of  the  test.  It  was  necessary  to  use  these  small  increments  since  the 
effects  of  creep  are  magnified  by  higher  temperature,  moisture,  or  stress 
levels.  In  these  cases,  the  initial  or  primary  portion  of  the  creep 
curve  is  changing  very  rapidly.  Thus,  it  is  necessary  to  keep  the  time 
increments  small  just  after  the  load  has  been  changed  or  incremented. 

The  time  increments  used  in  each  test  in  this  example  were  selected 
somewhat  arbitrarily.  It  is  clear  that  a  certain  amount  of  care  is  re¬ 
quired  when  selecting  the  time  increments.  While  the  60°C  test  required 
several  one  minute  increments,  the  -6.7°C  test  did  not.  Of  course  it 
is  always  possible  to  use  more  and  smaller  time  increments  for  any  par¬ 
ticular  loading  case.  However,  the  additional  accuracy  achieved  requires 
additional  computer  time  and  may  or  may  not  be  cost  effective.  Finer 
increments  than  those  used  in  Figure  12  were  tried,  but  no  significant 
improvement  in  accuracy  was  obtained. 

Uniaxial  creep-recovery  tests  were  also  simulated  for  the  Hercules 
epoxy  resin  data  generated  in  the  present  study.  Since  the  effects  of 
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temperature  and  the  parameters  necessary  to  describe  such  effects  were 
not  investigated,  two-hour  creep  followed  by  two-hour  recovery  tests  at 
various  stress  levels  were  modeled.  Results  for  2000  psi,  8000  psi,  and 
19,000  psi  uniaxial  tensile  simulations  at  21.0°C  and  zero  moisture  are 
shown  in  Figures  D-l  through  D-6  in  Appendix  D.  Figures  D-2,  D-4,  and 
D-6  show  the  creep  portion  of  the  curves  with  expanded  strain  scales. 

The  predictions  (symbols)  compare  well  with  the  experimentally  measured 
response  (solid  lines).  For  all  cases,  the  error  is  less  than  6  percent 
and  the  general  shape  of  the  predicted  response  matches  the  experimental 
response  quite  well.  Closer  examination  of  these  figures  reveals  that 
the  error  present  is  essentially  constant  for  all  times,  meaning  that  the 
error  was  generated  during  the  first  few  time  increments.  Again,  to 
keep  this  error  small,  as  shown  in  these  figures,  the  time  increments 
were  judiciously  chosen  for  the  primary  stages  of  creep  (when  the  load 
had  just  been  changed),  especially  at  higher  stress  levels.  As  expected, 
the  error  increases  with  higher  applied  stresses. 

While  the  recovery  predictions  for  all  simulations  (including  the 
Beckwith  data)  matched  the  expected  response  quite  well,  there  was  con¬ 
siderable  error  for  the  recovery  predictions  for  the  Hercules  resin  sub¬ 
ject  to  19,000  psi  (see  Figure  D-5).  No  reasonable  explanation  for  this 
anomaly  can  be  offered  at  this  time. 

Stress  contour  plots  for  the  above  simulation  have  not  been  presented 
since  the  resin  is  obviously  in  a  uniform  stress  state.  Element  stresses 
were  checked  immediately  after  application  of  a  load  and  at  various  times 
later  to  make  certain  that  there  were  no  changes.  Also,  the  element 
stresses  were  checked  to  make  certain  they  were  zero  upon  removal  of  the 
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load  and  at  various  times  following. 

Creep  of  a  Composite 

Using  the  data  for  the  Shell  epoxy  resin  generated  by  Beckwith  [12] 
and  the  properties  for  S2  glass  given  in  Table  I,  two-hour  creep  and  two- 
hour  recovery  tests  for  a  1000  psi  transverse  normal  stress  were  simulat¬ 
ed.  The  assumed  properties  of  the  fibers  were  based  upon  available  experi¬ 
mental  data  as  referenced.  Since  the  transverse  properties  of  the  fibers 
(i.e.,  Et,  V£t,  vtt,  at)  are  not  well-characterized,  their  values  are 
estimates  based  upon  the  existing  literature.  In  all  of  the  following 
composite  models,  the  fibers  are  assumed  to  behave  time-independently . 

That  is,  the  fibers  themselves  do  not  show  any  viscoelastic  behavior. 
Longitudinal  creep  data  generated  by  Beckwith  [12]  and  Irion  [24]  for 
unidirectional  glass  and  graphite  composites  support  this  assumption. 

The  results  of  each  simulation  are  compared  with  those  for  actual 
composite  specimens  [12]  subject  to  three  different  constant  temperatures 
of  -6. 7°  C  (20°F),  23.9°C  (75°F),  and  60°C  (140°F) .  They  are  plotted  in 
Figure  14.  The  moisture  content  was  taken  to  be  zero  for  all  simulations 
since  no  data  were  available  [12]  to  account  for  this  effect.  While 
Beckwith  [12]  used  S-901  glass  fibers  instead  of  S2  glass,  the  difference 
between  the  fiber  properties  should  not  affect  the  results  significantly, 
as  will  be  explained  shortly. 

The  agreement  for  low  temperatures  is  excellent;  however,  there  is 
considerable  error  at  the  high  temperature  of  60.0°C.  The  time  incre¬ 
ments  were  the  same  as  the  ones  used  previously  for  the  epoxy  resin 
modeling.  While  this  could  be  responsible  for  some  of  the  error  (see 
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TIME  (HOURS) 

Figure  14.  Plot  of  strain  versus  time  for  S2  glass  composite  (using  Shell  58-68R  epoxy  resin) 
subject  to  a  2.1  hour  creep  and  a  2.1  hour  recovery  test.  Predicted  values  are 
shown  by  the  symbols  while  the  solid  lines  represent  the  experimental  response  [12], 
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Figure  13),  it  certainly  is  not  the  sole  cause.  Since  there  was  reason¬ 
ably  good  agreement  for  the  neat  resin  at  the  high  temperature,  and 
since  transverse  normal  loading  is  a  matrix-dominated  response,  it 
would  be  reasonable  to  assume  that  the  results  presented  in  [  12]  for 
this  particular  temperature  are  questionable.  That  is,  during  a  trans¬ 
verse  loading,  the  fiber  properties  do  not  play  a  significant  role  in 
the  strain  response  of  the  composite.  While  the  fiber  volume  content 
does  affect  the  response,  as  will  be  demonstrated  later,  the  time-indpend- 
ent  modulus  and  transverse  strain  are  typically  characteristic  of  the 
type  of  matrix  material  used  (in  this  case  epoxy  resin).  Thus,  it  is 
logical  to  assume  that  the  transverse  time-dependent  (creep)  response  of 
a  composite  will  also  be  characteristic  of  the  type  of  matrix  material 
used.  Further  examination  of  the  value  generated  by  Beckwith  at  this 
particular  temperature  reveals  some  discrepancy  between  it  and  the  values 
for  other  fiber  angles.  The  data  in  Reference  [12]  show  the  time- independ¬ 
ent  portion  as  well  as  the  time-dependent  portion  of  the  creep  compliance 
to  be  gradually  increasing  as  the  angle  of  orientation  of  the  load 
(with  respect  to  the  fibers)  increases.  In  the  Beckwith  data  [12],  how¬ 
ever,  there  is  a  drop  in  the  time-independent  portion  and  a  drastic  in¬ 
crease  in  the  time-dependent  portion  of  the  creep  compliance  for  the 
angle  of  90°.  Thus,  it  is  likely  that  the  prediction  made  using  the 
present  analysis  represents  more  closely  the  true  material  response 
than  the  reduced  data  presented  in  Reference  [12]. 

The  following  comparisons  with  the  data  presented  by  Irion  [24] 
tends  to  support  this  statement.  Using  data  for  Hercules  3501-6  epoxy 
resin  determined  in  this  study,  and  the  properties  for  S2  glass  and 
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Hercules  AS  graphite  fibers  given  in  Table  I,  creep  tests  were  simulated. 
The  S2  glass  composites  were  subjected  to  -20  and  -23  ksi  (compressive) 
transverse  stresses  in  the  x-direction.  The  predictions  are  plotted 
along  with  actual  experimental  data  in  Figures  15  through  17.  Similarly, 
an  AS  graphite/epoxy  composite  subjected  to  a  -27  ksi  (compressive) 
transverse  stress  was  simulated  and  the  comparison  shown  in  Figure  18. 

In  all  four  figures  the  scale  has  been  expanded,  magnifying  the  error 
between  the  predictions  and  the  experiment.  In  Figures  15  through  18, 
predictions  are  shown  by  the  symbols  while  the  experimental  results  are 
represented  by  the  solid  lines.  All  simulations  were  made  assuming  the 
temperature  to  be  21°C  and  the  moisture  content  to  be  zero.  While  actual 
test  specimens  may  have  contained  some  moisture  due  to  the  relative  am¬ 
bient  humidity,  no  data  were  available  to  determine  this  influence  upon 
the  creep  parameters  used  in  the  analysis.  Therefore  its  effect  was 
neglected. 

For  the  comparisons  mentioned  above,  two  predictions  were  actually 
made,  one  including  a  cure  simulation  and  one  without  a  cure  simulation. 
The  predictions  with  a  cure,  represented  by  the*  symbols  in  Figures  15 
through  18,  involved  simulating  a  cooldown  from  the  177 °C  cure  tempera¬ 
ture  to  room  temperature  (21°C).  This  was  done  using  six  time-independ¬ 
ent  temperature  increments  to  achieve  210C,  and  was  followed  by  a  140- 
hour  relaxation  period.  The  relaxation  time  was  arbitrary  and  was  in¬ 
cluded  in  order  to  simulate  what  would  happen  to  the  stresses  in  a  com¬ 
posite  after  it  had  been  cured  and  allowed  to  "stand"  for  several  hours. 
It  was  discovered  that  longer  relaxation  times  tend  to  produce  lower 
and  more  uniform  stress  states.  This  will  be  demonstrated  later  using 
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Figure  16.  Plot  of  strain  versus  time  for  an  S2  glass  composite  (using  Hercules  3501-6  epoxy  resin) 
subject  to  a  14  hour  creep  test.  Predicted  values  are  shown  by  the  symbols  while  the 
solid  line  represents  the  experimental  data  [24]. 
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Figure  18.  Plot  of  strain  versus  time  for  an  AS  graphite  composite  (using  Hercules  3501-6  epoxy  resin) 
subject  to  a  14  hour  creep  test.  Predicted  values  are  shown  by  the  symbols  while  the 
solid  line  represents  the  experimental  data  [24]. 
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contour  plots  of  the  composite  model  during  a  cure  simulation. 

After  this  relaxation  period,  the  load  was  applied  and  the  model 
allowed  to  creep.  Creep  strains  plotted  in  Figures  15  through  18  for 
the  "with  cure"  simulation  were  "zeroed"  by  subtracting  the  thermal 
strain  prior  to  loading  from  the  subsequent  total  strains.  In  examin¬ 
ing  these  figures  it  will  be  noted  that  the  creep  strains  for  the  sim¬ 
ulation  which  involved  a  cure  cycle  are  higher  than  those  for  the  sim¬ 
ulation  which  did  not  involve  a  cure  cycle.  The  reason  is  that  the  in- 
inital  internal  stresses  in  an  uncured  composite  model  are  zero,  while 
the  internal  stresses  in  the  cured  composite  are  very  high.  These  in¬ 
itial  stresses  are  mainly  compressive  and  are  caused  by  the  difference 
in  the  coefficients  of  thermal  expansion  between  the  fiber  and  the  matrix. 
In  the  transverse  plane  (see  Figure  5),  the  fiber  contracts  during  the 
cooling  process,  while  the  matrix  contracts  even  more.  With  the  presence 
of  a  fiber,  however,  this  matrix  contraction  is  restricted,  producing 
compressive  residual  stresses  in  the  matrix.  Thus,  the  application  of  a 
compressive  stress  serves  to  raise  the  internal  stress  state,  which 
tends  to  make  the  composite  creep  more.  Actually,  the  cause  of  these 
residual  stresses  in  the  composite  is  more  complicated  than  just  described. 
The  relative  stiffnesses  and  Poisson's  ratios  of  both  the  fiber  and  the 
matrix  material  also  play  a  role  in  determining  these  stresses.  Adams 
and  Monib  [29]  have  described  the  influence  of  these  parameters  in  detail. 

Both  "with"  and  "without  cure"  predictions  were  made  in  order  to 
obtain  a  "range"  in  which  the  experimental  data  might  lie.  In  other 
words,  the  subsequent  response  of  a  composite  is  a  function  of  how  long 
it  has  been  since  the  specimen  was  fabricated.  In  theory,  the  internal 
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stresses  should  be  very  small  and  uniform  after  a  very  long  relaxation 
time.  It  is  expected  that  the  actual  experimental  data  will  lie  some¬ 
where  within  this  range,  and  the  predictions  shown  in  Figures  15  through 
18  seem  to  confirm  this.  In  Figures  16  and  17,  the  predictions  (for  the 
-23  ksi  stress)  are  over  and  under,  respectively,  the  actual  measured 
data  for  the  glass/epoxy  composite,  while  the  experimental  data  for  the 
-20  ksi  stress  (Figure  15)  fall  within  the  predicted  range.  In  figure 
18,  the  predictions  for  the  graphite/epoxy  composite  overestimate  the  ex¬ 
perimentally  measured  response.  There  is  considerable  scatter  in  the  ex¬ 
perimental  data  found  in  Reference  [24].  The  errors  between  the  predic¬ 
tions  and  these  experiments,  as  shown  in  Figures  15  through  18,  are  of 
the  same  magnitude  as  the  experimental  scatter.  Thus,  the  agreement  be¬ 
tween  the  predictions  and  experiment  is  quite  reasonable  considering  the 
limited  abount  of  data  available.  It  is  important  to  note  that  in  each 
of  these  figures  the  strain  scale  has  been  enlarged,  magnifying  the 
errors  between  the  predictions  and  the  experimental  data.  By  picturing 
these  plots  with  a  full  scale,  the  relative  error  is  put  into  proper  per¬ 
spective. 

While  Irion  [24]  used  a  slightly  lower  cure  temperature  for  a  longer 
period  of  time,  it  is  expected  that  the  results  using  his  cure  would  be 
very  close  to  those  actually  modeled.  In  fact,  the  lower  cure  tempera¬ 
ture  would  tend  to  induce  lower  initial  residual  stresses.  Thus,  such 
a  cure  simulation  would  fall  within  the  range  already  predicted. 

Contour  plots  for  the  glass/epoxy  composite  for  both  the  "with"  and 
"without  cure"  cycles  are  shown  in  Appendix  D.  Figures  D-7a  and  D-7b 
show  octahedral  shear  stress  and  strain  contours  immediately  after  cool- 


56 


down  of  the  composite,  while  Figures  D-7c  through  D-7h  show  the  composite 
after  0.5,  2.1,  and  140  hours  of  relaxation.  Figures  D-7i  and  D-7j 
show  octahedral  shear  stress  and  strain  contours,  respectively,  for  the 
composite  just  after  application  of  the  load,  and  Figures  D-7k  and  D-7£ 
show  these  two  plots  for  the  composite  after  the  load  has  been  applied 
for  140  hours. 

It  is  important  to  note  that  the  octahedral  shear  stress  plots  for 
all  of  the  following  examples  have  been  "normalized".  That  is,  all  of 
the  contour  values  shown  in  the  plots  have  been  divided  by  the  yield 
value  of  the  octahedral  stress.  The  yield  values  for  the  epoxy  resin 
were  previously  determined  in  Reference  [27].  Thus,  a  contour  value 
greater  than  one  signifies  that  the  region  has  yielded. 

It  is  clear  from  Figures  D-7a,  D-7c,  D-7e,  and  D-7g  that  the  inter¬ 
nal  stresses  are  relaxing  from  the  cure  cycle,  as  expected.  While  the 
shapes  of  the  contour  lines  appear  to  be  the  same,  the  magnitude  of  each 
contour  line  has  been  reduced  with  time.  It  is  interesting  to  compare 
the  contour  values  of  the  octahedral  shear  stresses  in  Figures  D-7i  and 
D-7k.  After  140  hours  of  elapsed  time  after  application  of  the  load, 
the  contours  are  seen  to  have  shifted,  relieving  the  high  stress  state 
in  the  lower  right-hand  corner.  The  concentration  of  stresses  in  this 
corner  was  induced  by  loading  the  composite  in  the  x-direction  (see 
Figure  4b) . 

The  octahedral  shear  strain  plots  (Figures  D-7b,  D-7d,  D-7f,  and 
D-7h)  reveal  only  minute  changes  of  the  internal  strain  during  relaxation, 
as  exhibited  by  both  the  shapes  and  magnitudes  of  the  contour  lines. 

This  is  reasonable  since  a  change  in  temperature  creates  a  symmetric 
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loading  and  the  net  force  on  the  boundary  due  to  this  type  of  loading  is 
zero.  Figures  D-7j  and  D-7 l  illustrate  the  change  in  octahedral  shear 
strain  immediately  after  the  load  has  been  applied,  and  140  hours  later. 
During  this  time  the  strain  in  the  lower  right-hand  corner  is  seen  to 
substantially  Increase,  while  only  moderate  increases  in  strain  are  ob¬ 
served  for  the  upper  portions  of  the  model.  These  changes  in  strain  are 
responsible  for  the  redistribution  of  stresses  during  this  time  period. 

In  Figure  D-8,  plots  of  the  normal  stress  and  tangential  shear 
stress  on  the  fiber-matrix  interface  are  shown  for  the  above  simulation, 
immediately  after  cure,  140  hours  later,  immediately  after  application 
of  the  load,  and  140  hours  after  that.  Due  to  the  required  continuity 
of  the  normal  and  tangential  shear  stress  at  the  interface,  these  figures 
represent  the  interface  stress  state  in  both  fiber  and  matrix.  Since  a 
square  fiber  packing  was  assumed  in  this  model,  the  interface  stress 
plots  must  be  symmetric.  The  normal  stress  (Figures  D-8a  and  D-8c)  is 
symmetric  about  the  horizontal,  vertical,  and  45°  axes,  while  the  tan¬ 
gential  shear  stress  (Figures  D-8b  and  D-8d)  is  symmetric  about  the 
22*5°  and  67^°  axes.  The  shear  stress  is  antisymmetric  about  the  45° 
axis,  although  the  polar  plots  tend  to  distort  this  fact.  These  axes  of 
symmetry  will  change,  however,  for  a  nonsymmetrical  loading  (see  Figures 
D-8e  through  D-8h).  Comparing  Figures  D-8a  and  D-8b  with  D-8c  and  D-8d, 
and  D-8e  and  D-8f  with  D-8g  and  D-8h,  it  is  seen  that  there  is  a  reduc¬ 
tion  in  both  the  normal  and  tangential  shear  stresses  with  time.  For 
the  symmetrical  loading  caused  by  the  cooldown,  the  contour  lines  do  not 
change  in  shape  but  only  in  magnitude,  as  expected.  On  the  other  hand, 
there  is  a  slight  change  in  shape  as  well  as  magnitude  with  time  for  the 
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loaded  composite.  These  changes  are  in  agreement  with  the  upward  redis¬ 
tribution  of  stress  observed  in  Figures  D-7i  and  D-7k. 

Figure  D-9  represents  the  glass/epoxy  composite  in  which  curing  was 
not  simulated.  The  first  two  plots  are  octahedral  shear  stress  and  octa¬ 
hedral  shear  strain  contour  plots  for  the  instant  immediately  after  appli¬ 
cation  of  the  load,  while  the  last  two  are  for  140  hours  after  application 
of  the  load.  Since  the  internal  stresses  are  lower  in  this  model  than  in 
the  previous  model  (see  Figures  D-7i  and  D-7k) ,  there  is  less  redistri¬ 
bution  of  the  stresses,  as  evidenced  by  Figures  D-9a  and  D-9c. 

The  octahedral  shear  strain  (Figures  D-9b  and  D-9d)  is  again  seen 
to  increase  significantly  in  the  lower  right-hand  corner  (where  the  stress 
is  high),  with  only  moderate  increases  in  the  upper  portion  of  the  model. 
This  is  also  caused  by  the  lower  initial  stress  state  present  in  this 
composite  than  in  the  composite  with  the  cure.  Thus,  a  higher  stress 
state  in  a  composite  corresponds  to  a  greater  amount  of  creep  strain, 
which  in  turn  causes  a  greater  redistribution  of  the  stresses. 

The  normal  and  tangential  shear  stresses  on  the  fiber-matrix  inter¬ 
face  were  similar  to  those  shown  in  Figures  D-8e  through  D-8f.  While 
the  values  were  not  as  large  in  this  case,  the  changes  in  shape  and 
magnitude  were  very  similar. 

Effect  of  Fiber  Volume 

To  demonstrate  the  effect  of  fiber  volume,  two  identical  computer 
runs  were  made,  using  50  percent  and  63  percent  fiber  volume.  For  this 
example,  a  glass/epoxy  composite  was  subjected  to  a  -20  ksi  stress  at 
21°C  and  zero  moisture  for  140  hours.  The  composite' strain  response  is 
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plotted  in  Figure  19  for  the  two  fiber  volumes.  As  expected,  the  lower 
fiber  volume  content  composite  showed  more  creep  strain  than  did  the 
higher  fiber  volume  composite.  This  is  due  primarily  to  the  fact  that 
there  is  more  matrix  material  in  the  50  percent  fiber  volume  composite. 
This  is  consistent  with  the  results  presented  for  the  neat  resin  model 
(see  Figure  D-6) ,  which  showed  even  more  creep  than  did  the  50  percent 
fiber  volume  composite. 

Contour  plots  of  the  50  percent  fiber  volume  composite  are  shown  in 
Figure  D-10  for  times  immediately  after  application  of  the  load  and  140 
hours  later.  Comparison  of  these  plots  with  those  for  the  63  percent 
fiber  volume  composite  (Figure  D-9)  reveals  the  higher  octahedral  stress 
state  just  after  load  application  and  140  hours  later. 

The  octahedral  shear  strain  plots  for  the  50  percenc  fiber  volume 
composite  show  a  significant  change  in  the  upper  portion  of  the  model, 
just  as  the  higher  fiber  volume  composite  does.  This  indicates  that  a 
corresponding  rearrangement  of  the  stresses  is  occurring  in  this  area. 
Figures  D-lOa  and  D-lOc  illustrate  this  fact. 

The  normal  and  tangential  shear  stresses  along  the  fiber-matrix 
interface  (Figure  D-ll)  indicate  a  change  in  shape  similar  to  that  of 
the  63  percent  fiber  volume  composite.  The  magnitude  of  this  change 
is  much  smaller  for  the  50  percent  fiber  volume  composite,  however. 

This  means  that  while  there  is  a  greater  redistribution  of  stresses  in 
this  model,  they  are  occurring  away  from  the  fiber-matrix  interface. 

Relaxation  After  Cure 

Two  runs  were  made  in  order  to  examine  the  change  in  internal  stress 
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state  after  curing  the  composite.  Figures  D-7a  through  D-7h  demonstrate 
a  symmetric  distribution  of  the  contour  lines,  and  show  that  there  are 
no  changes  in  shape  of  these  lines  over  a  period  of  time.  This  is  rea¬ 
sonable  since  the  temperature,  assumed  to  be  uniform  throughout  the  com¬ 
posite,  would  create  a  symmetrical  stress  state  after  the  curing  process. 
While  the  shape  of  the  contour  plots  remained  unchanged ,  the  magnitudes 
of  the  lines  decreased  in  the  octahedral  shear  stress  plots.  This  indi¬ 
cates  that  the  relaxation  continues  to  occur  long  after  the  cooldown  pro¬ 
cess.  The  contour  plots  illustrate  that  the  majority  of  relaxation  takes 
place  within  a  few  hours,  however. 

The  results  of  a  modified  cure  cycle  are  presented  in  Figure  D-12. 

A  composite  was  cooled  from  177°C  to  99°C,  held  for  2.1  hours,  then 
cooled  to  room  temperature  (21°C)  and  allowed  to  relax  for  another  2.1 
hours.  Figures  D-12a  and  D-12b  show  contours  for  the  time  immediately 
following  the  cooldown  to  99°C,  while  Figures  D-12c  and  D-12d  illustrate 
the  stress  state  after  2.1  hours  of  relaxation.  In  Figures  D-12e  and 
D-12f,  the  contour  plots  are  shown  immediately  after  the  cooldown  to 
21°C,  and  Figures  D-12g  and  D-12h  for  2.1  hours  later. 

A  slight  reduction  in  the  octahedral  shear  stresses  was  observed 
in  this  case  relative  to  those  for  the  composite  cured  without  the 
intermediate  hold  (see  Figures  D-7e  and  D-7f).  This  indicates  that  it 
is  possible  to  "slow  cool"  a  composite  and  achieve  a  lower  internal 
stress  state  at  the  end  of  the  cooldown.  By  changing  the  intermediate 
temperature  and/or  the  length  of  time  the  composite  is  held  at  this 
temperature,  an  "ideal"  curing  cycle  may  be  found  using  the  current 
analysis.  This  can  be  extremely  beneficial  since  it  is  usually 


62 


desirable  to  have  a  composite  in  the  lowest  stress  state  possible  follow¬ 
ing  a  cure.  This  analysis  could  also  be  used  to  simulate  the  "annealing*' 
of  a  unidirectional  laminate.  In  other  words,  if  the  composite  were  re¬ 
heated  following  a  cooldown  and  slowly  cooled  again,  this  analysis  could 
be  used  to  determine  a  suitable  cooling  profile  to  reduce  the  residual 
stresses  in  the  composite. 

The  octahedral  shear  strains  in  this  example  case  did  not  change 
significantly  with  time,  due  to  the  symmetry  and  type  of  loading.  As 
mentioned  previously,  a  thermal  loading  produces  a  symmetric  Stress  state 
with  no  net  stresses  on  the  boundary.  Therefore,  the  total  strain  ex¬ 
hibited  by  the  composite  will  be  very  small.  The  interface  stresses 
changed  very  little  throughout  this  curing  cycle,  and  thus  are  not 
presented. 

Cyclic  Loading 

To  simulate  a  cyclic  loading,  a  glass/epoxy  composite  was  subjected 
to  a  plus  and  minus  5000  psi  uniaxial,  x-direction  loading.  No  cure  was 
simulated  so  that  the  effects  of  the  cyclic  loading  could  be  isolated  and 
studied  independently.  The  stress  was  applied  for  2.1  hours,  then  re¬ 
versed  for  another  2.1  hours.  About  one  and  one-half  cycles  were  modeled, 
and  the  total  strain  results  are  plotted  in  Figure  20.  Here  it  is  easily 
seen  that  the  maximum  absolute  value  of  the  strain  just  prior  to  the  re¬ 
versal  of  the  stress  is  constant  for  each  half  cycle.  At  higher  values 
of  stress,  it  is  expected  that  this  value  would  not  remain  constant.  To 
illustrate  this  point,  the  same  composite  was  subjected  to  a  plus  and 
minus  20  ksi  x-direction  loading,  a  plot  of  strain  versus  time  being 
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shown  in  Figure  21.  Here  the  maximum  absolute  value  of  the  strain  just 
after  each  stress  reversal  actually  increases  with  time.  That  is,  the 
amount  of  time-independent  strain  produced  by  the  first  reversal  of  _ 
stress  at  2.1  hours  is  greater  than  twice  the  amount  of  strain  produced 
by  the  application  of  the  initial  20,000  psi  stress  at  zero  hours.  This 
is  caused  by  the  program  using  a  different  stiffness  for  the  matrix  ma¬ 
terial  during  the  stress  reversal.  It  would  be  expected  that  an  incre¬ 
ment  of  -40,000  psi  would  produce  twice  the  amount  of  strain  that  the 
+20,000  psi  stress  did  if  the  material  behavior  were  elastic.  However, 
the  program  which  assumes  the  material  behavior  to  be  linear  initially, 
now  considers  the  matrix  to  be  yielded.  Since  it  cannot  anticipate  the 
next  loading  increment,  it  assumes  the  material  to  behave  in  a  nonlinear 
(plastic)  manner  and  hence  calculates  the  corresponding  tangent  modulus. 
This  modulus  will  be  lower  than  the  elastic  modulus,  which  gives  rise  to 
the  greater  amount  of  strain.  This  problem  could  have  been  avoided  by 
inputing  a  smaller  increment  of  stress  (-1000  psi  for  example)  followed 
by  a  larger  increment  of  stress  (-39,000  psi).  By  doing  this,  the  pro¬ 
gram  would  have  unloaded  the  composite  plastically  during  the  first 
increment,  sensed  that  the  composite  was  being  unloaded,  and  unloaded  the 
composite  elastically  during  the  second  increment.  Actually,  to  accurate¬ 
ly  model  stress  reversal  at  this  stress  level,  several  increments  are 
needed  to  apply  and  reverse  the  20,000  psi  stress  during  the  time-inde¬ 
pendent  increments.  This  was  not  done  in  this  case  since  there  was  only 
a  limited  amount  of  epoxy  matrix  creep  data  available  anyway.  Rather, 
the  idea  here  was  to  demonstrate  the  capability  of  the  analysis. 

It  is  also  interesting  to  note  that  the  strain  due  solely  to  creep 


Plot  of  strain  versus  time  for  a  glass/epoxy  composite  subject  to 
a  +20,000  psi  stress  for  2.1  hour  intervals. 


actually  Increases  with  each  half  cycle.  This  is  caused  by  the  fact 
that  there  are  two  pseudo  applied  stresses  immediately  after  a  stress 
change  which  contribute  to  the  total  creep  strain.  Since  these  two 
stresses  have  just  been  applied,  their  (primary  creep)  contributions  are 
much  larger  than  those  for  stresses  which  have  been  applied  for  some 
time  (secondary  creep).  As  an  example,  for  time  2.1<t<4.2  hours,  three 
contributions  of  strain  are  present  (see  Appendix  B) .  These  contributions 
are  from  the  +20,000  psi  stress  applied  at  time  t =0.0  hours  and  two  -20,000 
psi  stresses  applied  at  time  t=2.1  hours.  Examination  of  the  data  reveals 
that  the  creep  strains  due  to  the  pseudo  application  of  the  two  -20,000 
psi  stresses  at  2.1  hours  are  much  greater  than  those  due  to  the  +20,000 
psi  stress  applied  at  the  beginning  of  the  test.  The  contributions  of 
strain  due  to  the  +20,000  psi  stress  are  secondary  creep  effects  while 
those  for  recovery  and  the  application  of  the  -20,000  psi  stress  are 
primary  creep  effects  (see  Figure  lc).  Thus,  the  composite  shows  more 
creep  strain  during  the  second  half  cycle  than  it  did  during  the  first. 
Likewise,  the  pseudo  application  of  the  two  +20,000  psi  stresses  at  4.2 
hours,  plus  the  +20,000  psi  stress  applied  at  0.0  hours,  tends  to  negate 
the  effects  of  the  two  -20,0  psi  stresses  applied  at  2.1  hours. 

Contour  plots  of  the  internal  stress  state  immediately  after  appli¬ 
cation  and  just  prior  to  the  reversal  of  the  stress  for  each  half  cycle 
are  shown  in  Figure  D-13.  The  contour  plots,  as  before,  indicate  a  lower¬ 
ing  of  the  octahedral  shear  stress  with  time  after  the  load  has  just  been 
Inr resented.  While  both  the  octahedral  shear  stress  and  octahedral  shear 


•train  contour  shapes  changed  very  little,  the  magnitudes  do  indicate  a 
» v* line  of  the  stresses,  times  2.1  hours  later  being  shown.  It  is 
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interesting  to  observe  the  difference  in  contour  shapes  for  the  two  posi¬ 
tive  stress  cycles  during  the  test  (Figures  D-13a  through  D-13d  and  D-13i 
through  D-13£).  Immediately  after  application  of  the  second  +20,000  psi 
stress,  the  contour  values  are  higher  than  they  were  after  application  of 
the  initial  +20,000  psi  stress.  However,  two  hours  after  the  second 
+20,000  psi  stress  is  applied,  the  stresses  rearrange  themselves  in 
Figure  D-13k  very  much  the  same  as  they  did  two  hours  after  application 
of  the  first  +20,000  psi  stress  (Figure  D-13c) .  In  each  case,  the  highly 
stressed  region  has  moved  upward  (in  the  y-direction)  from  the  lower  right- 
hand  corner.  Thus,  even  after  a  stress  reversal,  the  highly  stressed 
regions  tend  to  redistribute  themselves  in  order  to  smooth  or  even  out  the 
stress  gradients  present.  The  octahedral  shear  strain  plots  show  this  re¬ 
distribution  by  the  large  values  of  strain  present  in  these  areas. 

In  going  from  the  +20,000  psi  (Figure  D-13c)  to  the  -20,000  psi 
(Figure  D-13e)  value  of  stress  (during  the  time-independent  increment), 
the  contour  shapes  also  rearrange  themselves  so  that  the  high  stress  con¬ 
centration  is  moved  upward.  On  the  other  hand,  in  going  from  the  negative 
(Figure  D-13g)  to  the  posititve  (Figure  D-13i)  value  of  stress,  the  high 
stress  concentration  area  is  pulled  downward. 

Contour  plots  for  the  plus  and  minus  5000  psi  loading  did  not  change 
much  since  the  material  behavior  did  not  exceed  the  elastic  limit.  The 
contour  shapes  remained  essentially  the  same  and  the  only  real  indications 
of  creep  were  the  changes  in  the  octahedral  shear  strain  plots. 

Plots  of  the  normal  and  tangential  shear  interface  stresses  for  the 
20  ksi  cycling  case  are  shown  in  Figure  D-14.  Again,  each  plot  corres¬ 
ponds  to  the  stress  state  just  before  and  after  the  load  is  incremented. 
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While  the  normal  interface  stresses  are  reduced  with  time  during  the 
tensile  loads,  they  increase  with  time  during  the  compressive  loading. 

On  the  other  hand,  the  tangential  shear  stresses  behave  in  the  opposite 
manner.  It  will  also  be  noted  that  the  signs  of  both  interface  stresses 
change  when  the  load  is  reversed. 

Multiaxial  Loading 

A  5000  psi  stress  was  applied  in  both  the  x-  and  y-directions  and 
the  composite  was  allowed  to  creep  for  140  hours.  Again,  no  cure  was 
simulated  so  that  the  effects  of  the  multiaxial  loading  could  be  isolated 
and  studied.  A  plot  of  x-direction  strain  versus  time  is  shown  in  Figure 
22.  The  y-direction  strain  versus  time  plot  is  not  presented  since  it  is 
identical  to  that  for  the  x-direction.  There  is  some  irregularity  In 
Figure  22  due  to  the  limited  number  of  points  plotted.  In  this  computer 
plot,  the  data  points  were  connected  with  straight  lines,  giving  the  false 
illusion  of  a  depression  between  14  and  45  hours. 

Contour  plots  of  the  internal  stress  state  are  shown  in  Figure  D-15. 

As  expected  for  a  symmetrical  loading,  little  redistribution  of  the  octa¬ 
hedral  shear  stresses  or  strains  occurred  even  after  140  hours.  In  fact, 
there  were  essentially  no  changes  in  the  magnitudes  of  the  octahedral 
shear  stress  at  all.  The  plots  of  the  octahedral  shear  strain  did  indicate 
that  creep  in  the  form  of  an  overall  expansion  of  the  model  was  taking 
place,  however. 

A  slightly  different  multiaxial  loading  case  was  also  simulated. 

Here,  a  15,000  psi  and  a  5000  psi  stress  were  applied  to  the  same  model 
in  the  x-  and  y-directions,  respectively,  for  140  hours.  Strain  versus 
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Figure  22.  Plot  of  strain  versus  time  (x-direction)  for  a  glass/epoxy  composite  subjected 

to  a  5000  psi  stress  in  the  x-direction  and  a  5000  psi  stress  in  the  y-direction 
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Figure  23.  Plot  of  strain  versus  time  (x-direction)  for  a  glass/epoxy  composite  subjected  to 
15,000  psi  stress  in  the  x-direction  and  a  5000  psi  stress  in  the  y-direction. 
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Figure  24.  Plot  of  strain  versus  time  (y-direction)  for  a  glass/epoxy  composite  subjected  to 
a  15,000  psi  stress  in  the  x-direction  and  a  5000  psi  stress  in  the  y-direction. 
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time  plots  are  shown  for  the  x-  and  y-dlrections  in  Figures  23  and  24. 

Note  that  while  the  x-direction  displayed  a  reasonable  amount  of  creep 
for  the  applied  stress,  the  y-direction  showed  very  little.  It  can  be 
seen  that  this  composite  exhibited  10  to  15  times  more  creep  in  the 
x-direction  than  in  the  y-direction.  By  comparing  the  strain  versus 
time  plots  for  the  direction  in  which  the  5000  psi  stress  was  applied  in 
both  multiaxial  loading  cases  (Figures  22  and  24) ,  it  is  seen  that  the 
specimen  with  the  symmetric  loading  exhibited  3  to  4  times  the  amount  of 
creep  strain  compared  to  the  specimen  which  was  not  symmetrically  loaded. 
This  illustrates  the  time-dependent  effect  of  a  constant  Poisson's  ratio. 

If  this  ratio  were  assumed  to  change  with  time,  the  response  would  ob¬ 
viously  be  different. 

Contour  plots  of  octahedral  shear  stress  for  the  nonsymmetrical  load¬ 
ing  case  (Figure  D-16)  show  a  slight  reduction  in  the  internal  stress 
state  after  140  hours.  On  the  other  hand,  the  octahedral  shear  strain 
contour  values  are  seen  to  increase  with  time.  As  in  the  symmetrical 
loading  case,  both  the  octahedral  shear  stress  and  strain  contour  shapes 
remain  essentially  unchanged. 

Interface  normal  and  tangential  shear  stress  plots  for  the  symmetric 
and  nonsymmetric  loading  cases  are  presented  in  Figures  D-17  and  D-18, 
respectively.  While  no  changes  in  the  tangential  shear  stresses  with 
time  for  either  simulation  were  observed,  the  normal  stresses  on  the  inter¬ 
face  Increased  slightly  for  the  symmetric  loading  and  decreased  for  the 
nonsymmetric  loading.  The  reason  for  this  is  that  the  applied  stresses  are 
causing  the  matrix  material  to  creep,  while  it  is  assumed  that  the  fiber 
does  not.  For  the  syimnetrical  loading  case,  the  stresses  caused  by  the 
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creep  strains  in  the  matrix  cannot  redistribute  themselves  and  therefore 
increase.  This  situation  is  analogous  to  placing  a  piece  of  material 
between  two  rigid  supports  and  raising  the  temperature.  As  the  material 
attempts  to  expand,  the  rigid  supports  restrict  the  expansion,  inducing 
higher  and  higher  stresses  in  the  material.  In  the  case  of  the  non- 
symmetric  loading,  the  stresses  can  redistribute  themselves,  as  they  did 
in  the  uniaxial  loading  cases,  in  order  to  relieve  stress  concentrations 
or  highly  stressed  areas. 

In  general,  it  can  be  said  that  for  the  composite  models  subjected  to 
a  decrease  in  temperature,  the  octahedral  shear  strain  plots  showed  only 
small  changes  with  time,  while  the  octahedral  shear  stress  plots  illus¬ 
trated  symmetrical  decreases.  On  the  other  hand,  models  subjected  to 
symmetrical  external  loadings  showed  balanced  increases  in  octahedral 
shear  strain  with  time,  with  only  minute  changes  in  octahedral  shear  stress. 
In  all  symmetric  loading  cases,  only  the  magnitudes  and  not  the  shapes  of 
the  contours  changed  with  time,  as  expected.  For  the  nonsymmetric  external 
loading  cases,  the  contour  lines  for  the  octahedral  shear  stress  plots  were 
reduced  in  magnitude,  and  redistributed,  in  order  to  relieve  highly  stress¬ 
ed  regions.  At  the  same  time,  the  octahedral  shear  strain  plots  illus¬ 
trated  this  rearrangement  of  the  stresses,  by  marked  increases  in  strain 
in  these  areas. 

Experimental  data  were  not  available  for  comparison  for  either  the 
load  cycling  or  the  multiaxial  creep  tests.  While  the  analysis  is  capable 
of  predicting  the  composite  response  due  to  either  of  these  loading  his¬ 
tories,  experimental  verification  would  be  useful.  A  series  of  uniaxial 
load  cycling  tests  should  be  run  at  various  stress  levels  in  order  to 
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confirm  the  behavior  predicted  by  the  analysis.  It  would  be  informative 
to  compare  the  amount  of  creep  strain  increase  per  half  cycle  as  predict¬ 
ed  by  the  analysis  (Figures  20  and  21)  with  that  of  experiment.  Multi- 
axial  creep  tests,  while  being  more  difficult  to  run,  could  be  used  to 
confirm  the  amount  of  creep  strain  predicted  (see  Figures  22  through  24) 
by  the  analysis. 

An  assumption  was  made  in  the  beginning  of  this  study  concerning  the 
coupling  between  multiaxial  applied  stresses.  Specifically,  this  assump¬ 
tion  was  that  the  only  coupling  between  multiaxial  stresses  was  due 
entirely  to  Poisson  effects.  By  performing  multiaxial  creep  tests,  the 
validity  of  this  assumption  could  be  evaluated.  It  is  felt  that  this  was 
a  good  assumption  due  to  the  excellent  correlation  between  the  predictions 
and  experiment  obtained  in  the  uniaxial  loading  cases.  Although  the  com¬ 
posite  is  loaded  uniaxlally,  the  matrix  is  in  a  very  complicated  trlaxlal 
state  of  stress.  Thus,  the  multiaxial  stress  response  of  the  epoxy  is 
actually  being  tested. 

Verification  of  the  changes  in  internal  stress  states  for  various 
cooldown  cycles  would  also  be  helpful.  Such  changes  are  difficult  to 
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CHAPTER  VI 


SUMMARY  AND  CONCLUSIONS 


The  several  examples  presented  in  the  previous  chapter  demonstrate 
only  a  few  of  the  many  potential  applications  of  the  analysis.  It  would 
be  possible  to  combine  several  of  these  features,  such  as  multiaxlal 
cycling,  thermal  cycling,  etc.,  to  model  many  different  situations.  How- 
ever,  the  purpose  here  is  to  give  an  indication  of  the  possible  uses  and 
applications  of  the  program.  As  mentioned  previously,  experimental  veri¬ 
fication  of  these  capabilities  would  be  helpful.  Such  verification  would 
allow  one  to  use  the  finite  element  program  with  assurance  that  predic¬ 
tions  of  complicated  hygrothermal  and/or  loading  histories  were  correct. 

In  summary,  the  analysis  is  capable  of  modeling  any  unidirectional 
composite  subjected  to  longitudinal  and/or  transverse  normal  and/or  hygro¬ 
thermal  loadings.  A  finite  element  computer  program  is  used  to  model  the 
composite.  The  constituent  material  parameters,  which  are  functions  of 
temperature,  moisture,  and  time,  are  input  to  the  program  through  two 
subroutines.  Thus,  material  properties  can  readily  be  changed.  The 
finite  element  analysis  used  in  the  computer  program  is  based  upon  a  dis¬ 
placement  formulation,  constant  strain  triangular  elements  with  three 
degrees  of  freedom  being  used.  The  analysis  is  pseudo-three-dimensional 
in  nature,  utilizing  the  conditions  for  generalized  plane  strain.  Non¬ 
linear  (elastoplastic)  material  behavior  is  included  by  means  of  the 
tangent  modulus  method.  While  this  program  is  incremental  in  nature,  and 
requires  the  use  of  several  steps  to  load  the  composite  into  the  non- 

i 


76 


linear  range,  it  can  readily  model  matrix  materials  with  "flat"  stress- 
strain  curves,  such  as  aluminum.  For  example,  the  modeling  of  a  boron/ 
aluminum  composite  would  present  no  difficulty. 

Viscoelastic  behavior  is  incorporated  into  the  analysis  by  using 
nonlinear  constitutive  equations  developed  by  Schapery  [8].  Creep  strains 
are  treated  as  initial  strains  and  are  calculated  by  linearly  approxi¬ 
mating  the  value  for  the  creep  compliance  over  a  time  increment.  Thus, 
it  is  necessary  for  these  increments  to  be  small  when  the  value  of  the 
creep  compliance  is  changing  rapidly.  While  uniaxial  loads  produce  com¬ 
plicated  triaxial  stress  states  in  the  composite,  it  is  assumed  that  the 
only  coupling  between  these  multiaxial  stresses  in  the  material  is  due 
to  Poisson  effects.  It  is  possible  to  include  higher  order  coupling 
effects  between  multiaxial  stresses;  however,  this  complication  has  been 
neglected  in  the  current  analysis.  Also,  Poisson's  ratio  has  been  assum¬ 
ed  to  be  time-independent  or  at  most  a  very  weak  function  of  time.  The 
validity  of  this  assumption  has  been  demonstrated  by  several  investiga¬ 
tors  [12,16,17].  By  making  use  of  it,  the  constitutive  equations  are 
simplified  immensely  (see  Appendix  A). 

The  present  analysis  solves  problems  which  involve  changes  of  stress, 
temperature,  and  moisture  content  with  time  as  a  series  of  step  loadings. 
That  is,  no  increments  of  load,  temperature,  or  moisture  are  permitted 
during  a  time  increment.  However,  increments  of  load,  temperature,  and / 
or  moisture  may  be  applied  simultaneously  using  a  time-independent  incre¬ 
ment.  Thus,  a  ramp  loading  must  be  simulated,  for  example,  by  a  series 
of  step  loads.  This  could  become  inefficient  as  the  number  of  steps 
grew  large,  due  to  the  nature  of  the  program.  As  explained  before,  each 
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successive  step  load  requires  two  additional  contributions  of  strain  to 
be  solved  for.  This  is  due  to  the  nonlinear  time-dependence  which  the 
matrix  material  displays. 

In  Chapter  III,  the  conditions  for  linear  viscoelastic  behavior 
were  discussed  and  a  test  for  determining  whether  a  material  is  linear 
or  nonlinear  was  described.  While  most  materials  display  some  linear 
behavior  at  low  stresses  and  temperatures,  these  same  materials  are 
usually  quite  nonlinear  at  higher  values  of  stress  and  temperature.  The 
transition  values  of  stress  and  temperature  (and  for  polymer-matrix  com¬ 
posites,  moisture)  are  different  for  each  material  and  must  be  determined 
by  experiment. 

The  nonlinear  viscoelastic  parameters  for  the  matrix  material  which 
are  input  to  the  analysis  can  be  determined  by  means  of  a  series  of  creep- 
recovery  tests  at  different  stress  levels,  as  demonstrated  in  Chapter  IV. 
To  include  the  effects  of  temperature  and  moisture,  a  series  of  creep- 
recovery  tests  at  various  temperatures  and  moisture  contends  must  be  run. 
Reduction  of  the  experimental  data  can  be  achieved  by  means  of  simple 
graphical  techniques. 

In  all  of  the  examples  presented  in  Chapter  V,  the  failure  criterion 
(present  in  the  existing  microraechanics  computer  program)  was  ignored. 
Since  failure  is  an  important  criterion  in  many  designs,  this  feature 
could  prove  very  useful.  First  failure,  as  currently  predicted  by  the 
analysis,  does  not  always  mean  that  a  composite  has  failed,  however. 

Many  composites  continue  to  carry  load  long  after  the  matrix  has  exper¬ 
ienced  many  microcracks.  For  this  reason,  crack  propagation  should  be 
included  in  future  analyses.  The  basic  micromechanics  program  has  been 
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modified  by  Murphy  and  Adams  [30]  to  include  crack  propagation.  Thus, 
this  capability  could  easily  be  incorporated  into  the  present  analysis. 

Other  future  considerations  might  include  the  addition  of  shear 
loading  into  the  analysis.  Currently,  work  is  underway  to  include  this 
capability  [31]  into  a  separate  micromechanics  program.  By  including 
this  type  of  loading,  all  time-dependent  parameters  necessary  for  use 
in  a  laminate  analysis  may  be  determined. 

A  3-D  finite  element  laminate  analysis  has  recently  been  developed 
by  Monib  [32],  The  program  uses  advanced  storage  and  computing  techni¬ 
ques  which  make  it  ideal  for  use  on  most  any  computer  system.  While  it 
does  not  include  viscoelastic  behavior  at  this  time,  it  could  easily  be 
modified  to  handle  such  effects.  This  analysis  could  then  be  used  to 
predict  the  viscoelastic  response  of  any  laminate  desired,  for  ultimate 
use  in  the  design  of  a  specific  component. 

An  important  feature  of  the  analysis  and  associated  computer  pro¬ 
gram  developed  in  the  present  study  is  that  only  constituent  material 
properties  need  to  be  evaluated.  This  avoids  the  difficulty  of  having 
to  evaluate  the  time-dependent  parameters  of  all  fiber-matrix  combin¬ 
ations.  Usually  only  matrix  materials  need  to  be  tested  for  time-depend¬ 
ence  since  most  fibers  do  not  show  (or  show  very  little)  viscoelastic 
behavior.  Thus,  a  designer  can  model  the  time-dependent  behavior  of 
many  different  types  of  composite  systems  with  only  a  small  amount  of 
experimental  data.  This  was  demonstrated  in  Chapter  V  for  glass/epoxy 
and  graphite/epoxy  composites.  Here,  only  the  epoxy  had  to  be  tested 
In  order  to  determine  the  viscoelastic  response  of  two  different 
composites. 
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In  general,  the  results  obtained  for  uniaxial  loading  compare  very 
well  with  the  experimental  data  available.  Comparisons  with  the  data 
generated  by  Irion  [24]  proved  the  analysis  to  be  within  experimental 
error.  While  only  a  limited  amount  of  data  was  available  for  comparison, 
it  is  felt  that  future  comparisons  will  be  equally  accurate.  The  analy¬ 
sis  is  capable  of  modeling  almost  any  type  of  loading  history  desired 
since  all  stress  components  are  accounted  for  individually.  Also,  it 
has  been  shown  that  it  is  possible  to  determine  an  ideal  curing  cycle 
for  a  unidirectional  composite,  or  the  internal  stress  state  of  a  com¬ 
posite  after  any  curing  cycle,  using  the  present  analysis.  Any  type 
of  thermal  or  moisture  history  can  be  modeled  to  determine  the  effects 
of  environmental  exposure  on  a  composite.  Th.'s  can  be  useful  in  de¬ 
termining  whether  a  particular  component  has  suffered  any  degradation 
during  this  period.  Combinations  of  mechanical  and  environmental  load¬ 
ings  can  also  be  modeled  using  the  program  to  determine  if  the  environ¬ 
ment  intensifies  or  relieves  the  internal  stresses  in  the  composite. 

Thus,  it  is  clear  that  the  analysis  and  associated  computer  program  of 
the  present  study  are  ideal  for  determining  the  time-dependent  stress 
state  and  overall  response  of  any  unidirectional  composite  under  many 
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APPENDIX  A 


EVALUATION  OF  THE  CONSTITUTIVE  EQUATIONS 


Single  Step  Input 

For  an  Isotropic  material  under  Isothermal  conditions  and  loaded 
uniaxially,  the  constitutive  equation  is  given  by  Eq.  (4)  of  Chapter 
III  as 


e(t) 


goD0a(t)  +  gl 


ADCf-f')  dtg20(t:>l  dt' 
dt ' 


(4) 


where  reduced  time  is  defined  by  Eqs.  (5) 


Since  the  term  a0  is  a  function  of  stress  (and  possibly  of  temperature 

and  moisture  also) ,  and  since  this  parameter  is  required  to  be  constant 

throughout  the  time  increment,  it  may  be  brought  outside  of  the  integral. 

First  consider  the  term  jJjgf J  in  Eq.  (4).  If  a  body  is  sub- 

dx' 

jected  to  the  loading  shown  in  Figure  A-l,  the  stress  history  may  be 
written  in  terms  of  the  unit  step  function  [33]. 

Let 


o(t)  -  oaU(T’-a)  (A-l) 

where  t*  is  the  variable  of  integration,  a  is  the  time  at  which  the  load 
is  applied,  and  the  unit  step  function  is  defined  as  [33] 
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Figure  A-l.  Plot  of  stress  versus  time  for  a  single-step  input 


Remembering  that  g2  is  constant  over  time,  and  substituting  Eq.  (A-l), 
we  can  write 


d[g2o(t)]  _  ~  d[oaU(x'-a)] 

t)  - 


(A-3) 


dx '  dx ' 

Since  oa  is  also  constant  over  time,  and  recognizing  that  the  derivative 
of  the  unit  step  function  is  the  Dirac  delta  function  [33],  Eq.  (A-3) 
becomes 

d[g2o(t)] 


dr' 


g2oalS(T’-a) 


(A-4) 


V  = 


dt'  w  Jt 
ao  ®o 


time 

terms 

(Eqs. 

5)  we  get 
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r  .  (T 

'  dtl  „  tMT  , 
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®<J  ao  L  ao 
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0  (A-5) 

1  and 

(A-5) 

into 

Eq.  (4)  we 

get 

c(t)  -  g0D0oaU(xr-a)  +  g 


2oa6(x ’-a)dx' 


(A-6) 


Upon  Integration,  the  constitutive  equation  becomes 


e(t)  -  g0D0aa  +  81g2AD(£l£a)oa  (A-7) 

for  all  t>ta.  Note  that  the  time-dependent  portion  of  the  compliance  is 
a  function  of  the  ratio  of  (t-ta)  and  a0 ,  which  is  reduced  time. 

Multiple  Step  Inputs 

A  slightly  more  complicated  case  than  that  just  presented  is  the 
two-step  input  shown  in  Figure  A-2. 


Figure  A-2.  Plot  of  stress  versus  time  for  a  two-step  input 


In  order  to  write  a  constitutive  equation  for  t>t_,  we  must  consider  the 
loading  as  consisting  of  three  parts:  + oa  applied  at  time  t“0  and  re¬ 
maining  unchanged,  -oa  applied  at  time  t*ta  and  remaining  unchanged  there¬ 
after  throughout  the  remaining  lifetime,  and  +o^  applied  at  time  t«ta 
and  also  remaining  unchanged  thereafter  throughout  the  remaining  lifetime. 
The  corresponding  constitutive  equation  using  Eq.  (4)  is 


e(t)  =  g^D0(aa)  + 


gj  ^  AD(f-4")gJ 


'(OdCf'JdT*  + 


g^>0(-Oa)  +  _AD(¥-»,)gJ(-H 


a-JSCt'-a)^'  + 


(A-8) 


SnDn<°b)  +  g^ 


j  ^  AD(¥-¥* 


)g“(ab)«(T’-a)dx’ 


where  the  substitution  defined  by  Eq.  (A-4)  has  been  made.  Note  that 
the  superscripts  on  the  terms  gQ ,  g1 ,  and  ^  denote  which  stress  level 
they  are  associated  with.  Upon  integration  and  cancelling  terms,  Eq. 


(A-8)  becomes 


e(t)  =  gjADa-vpgJoa  -  gjAD('P-'F’ )g*oa  + 

gJVb  +  8^AD<’r-’P2)°b 


(A-9) 


where  the  reduced  time  ¥  is 


^  -  4  +- 

Jo  a°  )taa°  a°  4 


(A-10a) 
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ta  T*-t, 
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The  equations  (A-lOb)  and  (A-lOc)  have  been  evaluated  at  t'*0  and 
T'“ta  respectively,  corresponding  to  the  times  at  which  each  load  was 
applied.  This  is  because  the  delta  function  in  the  intergrand  in  Eq. 
(A-8)  is  zero  everywhere  except  at  values  of  t'  corresponding  to  the 
load  application. 

Substituting  Eqs.  (A-10)  into  Eq.  (A-9)  we  get 


e(t)  =  —  + 

i 


ta  ,  t-ta 


b  a, Jt-t 


ffa  ~  gjg2AD 


(A-ll) 


g0D0°b  +  8^82AD 


t-ta\ 
—  ab 
ao  1 


After  rearranging  terms,  Eq.  (A-ll)  becomes 


e(t)  - 

80®0°b 
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Iab  ‘ 

b  a  / 

g.S/Dj 

t-tj.1 

t-ta 
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ao 
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(A-12) 


Recovery 

If  Eq.  (A-12)  is  used  where  o^  is  equal  to  zero,  as  shown  by 
Figure  A-3,  then  the  constitutive  equation  for  the  recovery  portion  of 
the  creep  curve  can  be  written  as 


er<*>  -  g^2AD(^f  + 

\a„  a„ 


b  a 

"  g,g  AD 
1  2 


t-t. 


Ja* 


(A-13) 


When  stresses  are  zero  during  the  recovery  portion  of  the  creep  curve, 
the  constants  g^,  g2,  and  a£  are  all  equal  to  one,  indicating  that  the 


material  behavior  is  linear  viscoelastic.  With  these  substitutions. 
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time,  t 

Figure  A-3.  Plot  of  stress  versus  time  for  a  creep-recovery  input 


Eq.  (A-13)  becomes 

er(t)  =  gf[AD(^-+  t-ta)  -  AD(t-ta)]oa.  (A-14) 

ao 

While  it  is  possible  to  represent  the  transient  creep  compliance  AD() 
by  almost  any  function,  it  is  convenient  to  use  a  power  law  in  time, 

AD('P)  =  CV11  .  (A-15) 

Substitution  of  Eq.  (A-15)  into  Eq.  (A- 14)  gives 


Letting 


cr(t)  -  g*  C  [ef  +  t-ta)n  -  (t-ta)n]  . 
ao 


(A-16) 


(A-17) 


dropping  the  superscript  a,  substituting  Eq.  (A-17)  into  Eq.  (A-16),  and 

rearranging  we  get 


er(t)  -  Cgzt£  [ (l+aaX)n  -  (aaX)nlaa 


(A-18) 
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The  amount  of  transient  or  creep  strain  just  prior  to  the  removal  of  the 

load  at  time  ta  is  defined  as 

Aea  =  e(ta)  -  e0  *  gjg2C(ta)noa  (A-19) 

where  eQ  is  the  initial  or  time-independent  strain.  Substituting  Eq. 

(A-19)  into  Eq.  (A-18)  we  get 

er(t)  =  ^  [ (l+aaA)n  -  (aaX)n] .  (A-20) 

gl 
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APPENDIX  B 


COMPUTER  PROGRAM 


The  computer  program,  originally  developed  by  Miller  and  Adams  [14], 
uses  simplex  triangular  elements  and  is  based  upon  generalized  plane 
strain  conditions,  which  permit  a  pseudo-three-dimensional  analysis. 
Nonlinear  (elastoplastlc)  material  behavior,  modeled  by  Von  Mises  yield 
criterion  and  the  Prandtl-Reuss  flow  rule,  is  incorporated  by  means  of 
the  tangent  modulus  method.  This  technique  allows  for  the  modeling  of 
matrix  materials  which  have  "flat”  or  "nearly  flat"  stress-strain  curves, 
such  as  aluminum. 

Temperature-  and  moisture-induced  strains  are  included  as  are  ma¬ 
terial  properties  which  are  a  function  of  temperature  and  moisture.  The 
efficiency  of  the  current  micromechanics  computer  program  lies  within  a 
method  which  allows  all  external  tractions  to  be  applied  and  solved  for 
simultaneously.  This  method  also  permits  changes  in  temperature  and 
moisture  to  be  applied  simultaneously  with  changes  in  load  during  the 
same  increment.  This  avoids  the  problem  of  having  to  solve  for  each 
loading  component  independently  and  using  superposition  to  obtain  the 
total  solution. 

Time-dependence,  incorporated  by  means  of  initial  strains  [15], 
uses  the  constitutive  equations  presented  in  Chapter  III.  Application 
of  any  load,  temperature,  or  moisture  increments  must  occur  in  the  time- 
independent  steps.  That  is,  no  changes  in  load,  temperature,  or  moisture 
are  allowed  to  occur  during  a  time  increment. 
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A  block  diagram  of  the  computer  program  is  shown  in  Figure  B-l. 

When  this  diagram  is  compared  to  the  one  in  Reference  [14],  the  only 
significant  difference  is  the  additional  subroutine  which  contains  the 
viscoelastic  parameters.  This  feature  allows  the  user  to  merely  sub¬ 
stitute  a  modified  subroutine  in  order  to  simulate  a  different  visco¬ 
elastic  matrix  material. 

The  program  previously  required  the  modification  of  another  sub¬ 
routine  in  order  to  model  a  different  matrix  and/or  fiber  in  a  composite. 
With  the  addition  of  time-dependence,  there  are  now  two  subroutines 
which  require  modification  in  order  to  model  a  different  composite 
material. 

A  typical  computer  run  involves  initialization  of  several  parameters 
in  the  subroutine  MAIN.  Next,  mesh  data  which  includes  node  point 
coordinates,  element  numbers,  element  material  types,  and  the  node  boun¬ 
dary  conditions,  are  read  into  the  program  by  subroutine  GDATA.  The 
grid  used  in  the  program  is  set  up  for  40  percent  fiber  volume;  however, 
subroutine  MADJ  automatically  recalculates  the  node  point  coordinates  for 
any  fiber  volume  desired.  Subroutine  MAIN  then  calls  subroutine  LOAD, 
which  reads  the  first  increment  of  temperature,  moisture,  or  load.  This 
subroutine  in  turn  calls  FORMK,  which  assembles  the  global  stiffness 
matrix  using  the  elastic  material  properties  subroutine  ELPROP  if  the 
element  is  in  the  linear  range,  or  the  plastic  material  properties  sub¬ 
routine  PLPROP  if  the  element  is  in  the  nonlinear  (elastoplastic)  range. 
For  the  first  increment,  the  behavior  is  automatically  assumed  to  be 
linear.  Should  a  succeeding  increment  be  a  time  increment,  the  elastic 
material  properties  subroutine  ELPROP  is  called  (to  get  certain  time- 
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FORMK 

Calls  subroutines  which  calculate  element  stiffness 
matrices.  Assembles  the  element  stiffness  matrices 
into  global  form.  Also  applies  boundary  conditions 
read  in  GDATA. 
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Figure  B-l.  (continued)  Analysis  Procedure  Flowchart 
(continued  on  the  following  page) 


Figure  B-l .  (continued)  Analysis  Procedure  Flowchart 
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independent  parameters)  and  then  the  time-dependent  materials  property 
subroutine  VISCO  is  called.  This  subroutine  calculates  the  elemental 
stiffness  matrix  (Eq.  14)  and  the  incremental  creep  compliance. 

When  formulation  of  the  global  stiffness  matrix  is  complete,  sub¬ 
routine  RESL  calculates  incremental  node  point  forces  (Eq.  15)  for  tem¬ 
perature,  moisture  or  time  increments.  Inversion  of  the  global  stiffness 
matrix  and  solution  for  the  node  displacements  occurs  in  subroutine 
SOLVE.  The  incremental  stresses  are  then  calculated  in  subroutine 
STRESS.  This  subroutine  also  calculates  the  octahedral  shear  stress 
and  strain  and  checks  for  yielding  and/or  failure  of  an  element.  If  no 
elements  have  failed,  the  program  returns  to  MAIN  for  another  increment. 
This  process  is  repeated  until  all  increments  are  exhausted,  or  an  ele¬ 
ment  fails,  at  which  time  the  program  automatically  terminates. 

In  order  to  model  a  multistep  input  such  as  shown  in  Figure  6,  the 
history  of  each  load  step  must  be  kept  track  of  Independently.  A  scheme 
using  a  random  access  file  technique  was  developed  for  this  purpose.  By 
using  this  method,  large  arrays  necessary  for  storage  of  each  stress 
history  were  not  needed.  Thus,  it  was  only  necessary  to  increase  the 
core  requirements  of  the  previous  micromechanics  program  in  Reference 
[14]  by  some  two  percent. 

An  example  of  this  scheme  is  as  follows.  Initial  load,  temperature, 
and  moisture  changes  which  occur  in  several  time-independent  increments 
produce  certain  stresses  in  each  element.  As  these  increments  are 
applied,  the  elemental  stresses  are  accumulated  and  stored  in  a  file 
called  FILE  ONE.  When  the  first  time  increment  is  applied,  the  contents 
of  FILE  ONE  are  copied  into  another  file  called  FILE  FIFTY.  This  file 


is  saved  for  later  use  when  the  load,  temperature  or  moisture  is  again 
incremented.  As  succeeding  time  Increments  are  applied,  incremental 
stresses  due  to  creep  (if  any)  are  accumulated  in  FILE  ONE.  At  this 
point,  the  stresses  in  FILE  ONE  are  the  total  stresses  in  each  element. 
When  more  increments  of  load,  temperature  or  moisture  are  applied,  sim¬ 
ulating  a  step  input,  the  stresses  which  are  contained  in  FILE  FIFTY  are 
copied  into  another  file  called  FILE  TWO,  with  a  change  of  sign  on  each 
stress.  This  serves  to  simulate  the  pseudo  application  of  a  -oa  stress, 
as  discussed  relative  to  multistep  loadings  in  Chapter  III.  At  the 
same  time,  the  contents  of  FILE  FIFTY  are  copied  into  another  file  called 
FILE  THREE.  As  the  time- independent  increments  of  load,  temperature,  or 
moisture  are  applied,  the  stresses  are  accumulated  in  FILE  THREE.  When 
the  next  time  increment  is  applied,  the  contents  of  FILE  THREE  are  copied 
over  the  contents  of  FILE  FIFTY.  Thus,  for  each  element  an  initial  "+oa 
stress  history"  beginning  at  time  t=0  corresponding  to  FILE  ONE,  a  ”-oa 
stress  history"  beginning  at  time  t=ta  corresponding  to  FILE  TWO,  and  a 
"+Ojj  stress  history"  beginning  at  time  t=ta  corresponding  to  FILE  THREE 
is  recorded.  A  third  step  load  in  time  similar  to  that  shown  in  Figure 
6  at  say,  time  t*tj,  would  require  two  more  files,  FILE  FOUR  for  a  "-a^ 
stress  history"  occuring  at  time  t-t^,  and  FILE  FIVE  for  a  "4oc  stress 
history"  occuring  at  time  tMt^.  Similarly,  extra  step  loads  would  re¬ 
quire  the  addition  of  two  files  per  step. 

The  solution  for  a  time  increment  for  times  ta<t<tfc  in  Figure  6 
requires  incremental  creep  stresses  to  be  calculated  for  each  stress 
history;  in  this  case  there  are  three.  These  incremental  stresses  are 
accumulated  in  their  respective  files  and  all  files,  i.e.,  ONE,  TWO,  and 


THREE  are  summed  at  the  end  of  an  inrement  to  give  the  total  stress  for 
each  element.  This  procedure  Is  then  repeated  for  additional  Increments. 


APPENDIX  C 
EXPERIMENTAL  RESULTS 


The  following  Figures  (C-l  through  C-10)  are  actual  computer  plots 
of  creep-recovery  data  for  Hercules  3501-6  epoxy  resin.  All  plots  at 
each  stress  level  were  similar  and  for  brevity  only  two  of  the  three  plots 
for  each  stress  level  tested  are  shown.  While  the  data  are  relatively 
smooth  for  the  higher  stress  levels,  they  are  quite  Irregular  for  the 
lower  values  of  applied  stress.  This  is  due  primarily  to  the  lower  out¬ 
put  of  strain  shown  by  the  specimens  subjected  to  a  lower  stress.  A 
smaller  output  requires  the  extensometer  to  be  set  on  a  more  sensitive 
scale  and  consequently,  the  extensometer  is  also  more  sensitive  to 
voltage  fluctuations,  noise  levels,  etc.  The  plots  also  appear  quite 
irregular  due  to  the  limited  number  of  data  points  taken.  In  these  plots, 
a  straight  line  was  used  to  connect  the  data  points.  During  the  data 
reduction  described  in  Chapter  IV,  a  smooth  curve  should  be  drawn  between 
these  points. 

Figures  C-ll  through  C-15  are  plots  of  the  constitutive  equations 
using  the  parameters  gQ ,  gj ,  g2,  and  a<,  calculated  from  the  creep-recov¬ 
ery  data.  The  experimental  creep  data  sampled  by  the  computer  for  all 
three  tests  are  indicated  by  the  stars.  The  stress  level  given  in  the 
upper  right-hand  corner  is  the  average  of  the  applied  stress  for  all 
three  tests.  There  is  excellent  agreement  between  experiment  and  the 
constitutive  equations  except  at  the  15,000  psi  stress  level.  By 
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examining  Figure  10  in  Chapter  IV,  the  reason  for  this  discrepancy  is 
evident.  The  value  for  gQDQ  obtained  through  data  reduction  and  the 
value  predicted  by  the  curve  in  this  plot  are  not  the  same,  and  since  the 
constitutive  equation  plotted  here  uses  the  value  predicted  by  the  curve, 
there  will  obviously  be  some  disagreement.  The  reason  for  this  discrep¬ 
ancy  in  Figure  10  is  not  known  at  this  time;  however,  the  experimental 
data  are  suspect,  as  discussed  in  Chapter  IV. 
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Figure  C-l.  Experimental  creep-recovery  data  for  Hercules  3501-6  epoxy 
resin  subject  to  a  2000  psi  uniaxial  compressive  stress. 
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Figure  C-4. 
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Experimental  creep-recovery  data  for  Hercules  3501-6  epoxy 
resin  subject  to  a  6000  psl  uniaxial  compressive  stress. 
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Figure  C-5.  Experimental  creep-recovery  data  for  Hercules  3501-6  epoxy 
resin  subject  to  an  8000  psi  uniaxial  compressive  stress. 
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Figure  C-6.  Experimental  creep-recovery  data  for  Hercules  3501-6  epoxy 
resin  subject  to  an  8000  psi  uniaxial  compressive  stress. 
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Figure  C-7.  Experimental  creep-recovery  data  for  Hercules  3501-6  epoxy 
resin  subject  to  a  16,000  psi  uniaxial  compressive  stress. 
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Figure  C-ll.  Plot  of  strain  versus  time  for  Hercules  3501-6  epoxy  resin  subjected  to  a  one  hour 
creep  test.  Experimental  values  are  shown  by  the  stars  while  the  solid  line 
represents  the  constitutive  equation. 


Plot  of  strain  versus  time  for  Hercules  3501-6  epoxy  resin  subjected  to  a  one  hour 
creep  test.  Experimental  values  are  shown  by  the  stars  while  the  solid  line 
represents  the  constitutive  equation. 
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Figure  C-14.  Plot  of  strain  versus  time  for  Hercules  3501-6  epoxy  resin  subjected  to  a  one  hour 
creep  test.  Experimental  values  are  shown  by  the  stars  while  the  solid  line 
represents  the  constitutive  equation. 
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Figure  C-15.  Plot  of  strain  versus  time  for  Hercules  3501-6  epoxy  resin  subjected  to  a  one  hour 
creep  test.  Experimental  values  are  shown  by  the  stars  while  the  solid  line 
represents  the  constitutive  equation. 


APPENDIX  D 


NUMERICAL  RESULTS 


In  Figures  D-l  through  D-6,  predictions  made  for  the  Hercules  3501- 
6  epoxy  resin  using  the  computer  program  are  plotted  with  the  data  which 
were  input  to  the  program.  The  predictions  are  represented  by  the  sym¬ 
bols  (squares)  while  the  lines  show  the  actual  (input)  response.  The 
input  response  was  obtained  by  using  the  constitutive  equations  and  a 
hand  calculator.  Figures  D-2,  D-4,  and  D-6  are  replotted  from  Figures 
D-l,  D-3,  and  D-5,  respectively,  with  the  scales  enlarged.  This  helps 
to  magnify  any  errors  present,  which  incidently  are  less  than  6  percent 
in  all  cases. 

Figures  D-7  through  D-18  are  octahedral  shear  stress  and  strain, 
and  normal  and  tangential  interface  stress  plots  for  various  examples 
presented  in  Chapter  V.  It  should  be  mentioned  that  the  octahedral 
shear  stress  plots  have  been  "normalized"  by  dividing  the  contour  values 
by  the  octahedral  shear  yield  stress.  This  means  that  contour  values 
greater  than  one  indicate  that  the  region  has  yielded.  The  octahedral 
shear  yield  stress  is  a  function  of  temperature  and  moisture  and  was  de¬ 
termined  to  be  5397  psi  [27]  for  all  plots  shown  except  for  Figures 
D-12a  and  D-12c.  In  these  figures  the  temperature  is  99°C,  which  is 
different  from  all  other  models  where  the  temperature  is  21°C.  In  this 
case  the  octahedral  shear  stress  yield  value  is  2813  psi  [27]. 
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Expanded  scale  plot  of  strain  versus  time  for  Hercules  3501-6  epoxy  resin  subjected 
to  a  2.1  hour  creep  and  a  2.1  hour  recovery  test.  Predicted  values  are  shown  by 
the  symbols  while  the  line  represents  the  actual  experimental  response. 
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while  the  line  represents  the 


Figure  D-4.  Expanded  scale  plot  of  strain  versus  time  for  Hercules  3501-6  epoxy  resin  subjected 
to  a  2.1  hour  creep  and  a  2.1  hour  recovery  test.  Predicted  values  are  shown  by 
the  symbols  while  the  line  represents  the  actual  experimental  response. 
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Figure  D-6.  Expanded  scale  plot  of  strain  versus  time  for  Hercules  3501-6  epoxy  resin  subjected 
to  a  2.1  hour  creep  and  a  2.1  hour  recovery  test.  Predicted  values  are  shown  by 
the  symbols  while  the  line  represents  the  actual  experimental  response. 
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a)  Octahedral  shear  stress  b)  Octahedral  shear  strain 

elapsed  time:  0.0  hours  elapsed  time:  0.0  hours 

(immediately  after  cooldown  to  21°C) 


c)  Octahedral  shear  stress  d)  Octahedral  shear  strain 

elapsed  time:  0.5  hours  elapsed  time:  0.5  hours 


Figure  D-7.  Octahedral  shear  stress  (normalized)  and  octahedral  shear 
strain  plots  for  glass/epoxy  subjected  to  a  cooldown  from 
the  177*C  cure  temperature  to  21°C,  followed  by  a  140 
hour  relaxation  period,  and  then  a  -20  ksi  stress  for 
another  140  hours. 


g)  Octahedral  shear  stress  h)  Octahedral  shear  strain 

elapsed  time:  140  hours  elapsed  time:  140  hours 


Figure  D-7.  (continued)  Octahedral  shear  stress  (normalized)  and 

octahedral  shear  strain  plots  for  glass/epoxy  subjected 
to  a  cooldown  from  the  177°C  cure  temperature  to  21°C, 
followed  by  a  140  hour  relaxation  period,  and  then  a 
-20  ksi  stress  for  another  140  hours. 


Figure  D-8.  Interface  stresses  for  glass/epoxy  subjected  to  a  cool¬ 
down  from  the  177°C  cure  temperature  to  21*C  followed  by 
a  140  hour  relaxation  period,  and  then  a  -20  ksi  stress 
for  another  140  hours. 
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elapsed  tine:  140  hours  elapsed  time:  140  hours 

(immediately  after  application  of  o^  *  -20  ksi) 


elapsed  time:  280  hours 


elapsed  time:  280  hours 


Figure  D-8.  (continued)  Interface  stresses  for  glass/epoxy  subjected 
to  a  cooldown  from  the  177°C  cure  temperature  to  21 ®C 
followed  by  a  140  hour  relaxation  period,  and  then  a  -20 
ksi  stress  for  another  140  hours. 


a)  Octahedral  shear  stress  b)  Octahedral  shear  strain 

elapsed  time:  0.0  hours  elapsed  time:  0.0  hours 

(Immediately  after  application  of  ox  =  -20  ksi) 


c)  Octahedral  shear  stress 
elapsed  time:  140  hours 


d)  Octahedral  shear  strain 
elapsed  time:  140  hours 


Figure  D-9.  Octahedral  shear  stress  (normalized)  and  octahedral  shear 
strain  plots  for  glass/epoxy  subjected  to  a  -20  ksi  stresi 
for  140  hours  (no  prior  temperature  history). 


a)  Octahedral  shear  stress  b)  Octahedral  shear  strain 

elapsed  time:  0.0  hours  elapsed  time:  0.0  hours 

(Immediately  after  application  of  ax  -  -20  ksi) 


c)  Octahedral  shear  stress  d)  Octahedral  shear  strain 

elapsed  time:  140  hours  elapsed  time:  140  hours 


Figure  D-10.  Octahedral  shear  stress  (normalized)  and  octahedral 

shear  strain  plots  for  a  50  percent  fiber  volume  glass/ 
epoxy  composite  subjected  to  a  -20  ksi  stress  for  140 
hours  (no  prior  temperature  history) . 


a)  Normal  stresses  (ksi)  b)  Shear  stress  (ksi) 

elapsed  time:  0.0  hours  elapsed  time:  0.0  hours 

(Immediately  after  application  of  ov  *  -20  ksi) 


c)  Normal  stresses  (ksi) 
elapsed  time:  140  hours 


d)  Shear  stresses  (ksi) 

elapsed  time:  140  hours 


Figure  D-ll.  Interface  stresses  for  a  50  percent  fiber  volume  glass/ 
epoxy  composite  subjected  to  a  -20  ksi  stress  for  140 
hours  (no  prior  temperature  history) . 
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a)  Octahedral  shear  stress 
elapsed  time:  0.0  hours 


b)  Octahedral  shear  strain 
elapsed  time:  0.0  hours 


(immediately  after  cooldown  to  99°C) 


c)  Octahedral  shear  stress  d)  Octahedral  shear  strain 

elapsed  time:  2.1  hours  elapsed  time:  2.1  hours 


Figure  D-12.  Octahedral  shear  stress  (normalized)  and  octahedral  shear 
strain  plots  for  glass/epoxy  subjected  to  a  cooldown  from 
the  1778C  cure  temperature  to  99°C,  a  2.1  hours  relaxation 
period,  further  cool  from  99°C  to  21*C,  and  another  2.1 
hour  relaxation  period. 


e)  Octahedral  shear  stress  f)  Octahedral  shear  strain 

elapsed  time:  2.1  hours  elapsed  time:  2.1  hours 

(immediately  after  cooldown  to  21°C) 


g)  Octahedral  shear  stress  h)  Octahedral  shear  strain 

elapsed  time:  A. 2  hours  elapsed  time:  4.2  hours 


Figure  D-12.  (continued)  Octahedral  shear  stress  (normalized)  and 

octahedral  shear  strain  plots  for  glass/epoxy  subjected 
to  a  cooldown  from  the  177°C  cure  temperature  to  99°C, 
a  2.1  hour  relaxation  period,  further  cool  from  99eC  to 
21°C,  and  another  2.1  hour  relaxation  period. 
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a)  Octahedral  shear  stress 
elapsed  time:  0.0  hours 


b)  Octahedral  shear  strain 
elapsed  time:  0.0  hours 


(immediately  after  application  of  ax  =  +20  ksi) 


c)  Octahedral  shear  stress 
elapsed  time:  2.1  hours 


d)  Octahedral  shear  strain 
elapsed  time:  2.1  hours 


Figure  D-13.  Octahedral  shear  stress  (normalized)  and  octahedral 

shear  strain  plots  for  glass/epoxy  subjected  to  a  +20 
ksi  stress  for  2.1  hours,  a  -20  ksi  stress  for  2.1  hours, 
and  a  +20  ksi  stress  for  another  2.1  hours  (no  prior 
temperature  history) . 


e)  Octahedral  shear  stress 
elapsed  time:  2.1  hours 


f)  Octahedral  shear  strain 
elapsed  time:  2.1  hours 


(immediately  after  stress  reversal,  to  ox  =  -20  ksi) 
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i)  Octahedral  shear  stress  j)  Octahedral  shear  strain 

elapsed  time:  4.2  hours  elapsed  time:  4.2  hours 

(immediately  after  stress  reversal,  to  ox  =  +20  ksi) 


k)  Octahedral  shear  stress 
elapsed  time:  6.3  hours 


t)  Octahedral  shear  strain 
elapsed  time:  6.3  hours 


Figure  D-13.  (continued)  Octahedral  shear  stress  (normalized)  and 

octahedral  shear  strain  plots  for  Rlass/epoxv  subjected 
to  a  +20  ksi  stress  for  2.1  hours,  and  a  -20  ksi  stress 
for  2.1  hours,  and  a  +20  ksi  stress  for  another  2.1 
hours  (no  prior  temperature  history). 
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elapsed  time:  0.0  hours  elapsed  time:  0.0  hours 

(immediately  after  application  of  ox  =  +20  ksi) 


elapsed  time:  2.1  hours 


elapsed  time:  2.1  hours 


Figure  D-14.  Interface  stresses  for  glass/epoxy  subjected  to  a  +20 

ksi  stress  for  2.1  hours,  a  -20  ksi  stress  for  2.1  hours, 
and  a  +20  ksi  stress  for  another  2.1  hours  (no  prior 
temperature  history). 


-45  0 

e)  Normal  stresses  (ksl) 

elapsed  time:  2.1  hours 


-10  0 

f)  Shear  stresses  (ksl) 

elapsed  time:  2.1  hours 


(Immediately  after  stress  reversal,  to  =  -20  ksi) 


-45  0 

g)  Normal  stresses  (ksl) 
elapsed  time:  4.2  hours 


-K>  0 

h)  Shear  stresses  (ksi) 

elapsed  time:  4.2  hours 


Figure  D-14.  (continued)  Interface  stresses  for  glass/epoxy  subjected 
to  a  4-20  ksi  stress  for  2.1  hours,  a  -20  ksi  stress 
for  2.1  hours,  and  a  +20  ksi  stress  for  another  2.1  hours 
(no  prior  temperature  history) . 


a)  Octahedral  shear  stress  b)  Octahedral  shear  strain 

elapsed  time:  0.0  hours  elapsed  time:  0.0  hours 

(immediately  after  application  of  ux  =  ov  «  +  5  ksi) 


c)  Octahedral  shear  stress  d)  Octahedral  shear  strain 

elapsed  time:  140  hours  elapsed  time:  140  hours 


Figure  D-15.  Octahedral  shear  stress  (normalized)  and  octahedral 

shear  strain  plots  for  glass/epoxy  subjected  to  a  5  ksi 
stress  in  the  x-dlrection  and  a  S  ksi  stress  in  the  y- 
direction  (no  prior  temperature  history). 
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a) 


Octahedral  shear  stress 
elapsed  time:  0.0  hours 


b)  Octahedral  shear  strain 
elapsed  time:  0.0  hours 


(immediately  after  application  of  ox  =  +  15  ksi,  oy  =  +  5  ksi) 
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c)  Octahedral  shear  stress 
elapsed  time:  140  hours 


d)  Octahedral  shear  strain 
elapsed  time:  140  hours 


Figure  D-16. 


Octahedral  shear  stress  (normalized)  and  octahedral 
shear  strain  plots  for  glass/epoxy  subjected  to  a  15  ksi 
stress  in  the  x-direction  and  a  5  ksi  strsss  in  the  y- 
direction  (no  prior  temperature  history). 


